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FOREWORD 


This report is the fifth revision of the document Observation Model and Parameter ParUals for 
the JPL VLBJ Parameter Estimation Software * MODEST ’ - 1991, dated August 1, 1991, which it 
supersedes. A number of model revisions and improvements were made from 1991 to 1994. They 
are briefly enumerated in the abstract. The present document corresponds to MODEST version 278, 
which has been in use since June, 1994. The authors hope to publish revisions of this document m 
the future, as modeling improvements warrant. 
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ABSTRACT 


This report is a revision of the document Observation Model and Parameter Partials for the JPL 
VLBI Parameter Estimation Software * MODEST * — 1991 , dated August 1, 1991. It supersedes that 
document and its four previous versions (1983, 1985, 1986, and 1987). A number of aspects of the 
very long baseline interferometry (VLBI) model were improved from 1991 to 1994. Treatment of tidal 
effects is extended to model the effects of ocean tides on universal time and polar motion (UTPM), 
including a default model for nearly diurnal and semidiurnal ocean tidal UTPM variations, and partial 
derivatives for all (solid and ocean) tidal UTPM amplitudes. The time-honored U K\ correction” for 
solid Earth tides has been extended to include analogous frequency-dependent response of five tidal 
components. Partials of ocean loading amplitudes are now supplied. The Zhu-Mathews-Oceans- 
Anisotropy (ZMOA) 1990-2 and Kinoshita-Souchay models of nutation are now two of the modeling 
choices to replace the increasingly inadequate 1980 International Astronomical Union (IAU) nutation 
series, A rudimentary model of antenna thermal expansion is provided. Two more troposphere 
mapping functions have been added to the repertoire. Finally, correlations among VLBI observations 
via the model of Tteuhaft and Lanyi improve modeling of the dynamic troposphere. A number of 
minor misprints in Rev. 4 have been corrected. 
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SECTION 1 


INTRODUCTION 


In applications of radio interferometry to geodynamics and astrometry, observed values of group 
delay and phase delay rate obtained from observations of many different radio sources must be passed 
simultaneously through a multiparameter estimation routine to extract the significant model param- 
eters. As the accuracy of radio interferometry has improved, increasingly complete models for the 
delay and delay rate observables have been developed. This report describes the current status of the 
delay model used in the Jet Propulsion Laboratory multiparameter estimation program “MODEST”, 
which is the successor to the “MASTERFIT” code developed at JPL in the 1970s. It is assumed 
that the reader has at least a cursory knowledge of the principles of very long baseline interferome- 
try (VLBI). Some references which provide an introduction are the book by Thompson, Moran, and 
Swenson (1986), and two reports by Thomas (1981, 1987), 

The delay model is the sum of four major model components: geometry, clock, troposphere, 
and ionosphere. Sections 2 through 5 present our current models for these components, as well as 
their partial derivatives with respect to parameters that are to be adjusted by multiparameter fits 
to the data. The longest section (2) deals with the purely geometric portion of the delay and covers 
the topics of time definitions, tidal and source structure effects, coordinate frames, Earth orientation 
(universal time and polar motion), nutation, precession, Earth orbital motion, wave front curvature, 
gravitational bending, and antenna offsets. Section 6 describes the technique used to obtain the delay 
rate model from the delay model. Section 7 gives the values of physical constants used in MODEST, 
while Section 8 outlines model improvements that may be required by more accurate data in the 
future. 


1 



SECTION 2 

GEOMETRIC DELAY 


The geometric delay is that interferometer delay which would be measured by perfect instrumen- 
tation, perfectly synchronised, if there were a perfect vacuum between the observed extragalactic or 
Solar-System sources and the Earth-based instrumentation. For Earth-fixed baselines, this delay can 
be as large as 20 milliseconds, changing rapidly (by up to 3.1 per second) as the Earth rotates. In 
general the geometric component is by far the largest component of the observed delay. The main 
complexity of this portion of the model arises from the numerous coordinate transformations necessary 
to relate the reference frame used for locating the radio sources to the Earth-fixed reference frame in 
which station locations are represented. 

In the following we will assume, unless otherwise stated, that “celestial reference frame” means 
a reference frame in which there is no net proper motion of the extragalactic radio objects which 
are observed by the interferometer. This is only an approximation to some truly “inertial” frame. 
Currently, this celestial frame implies a geocentric, equatorial frame with the equator and equinox of 
2000 January 1.5 (J2000) as defined by the 1976 International Astronomical Union (IAU) conventions, 
including the 1980 nutation series (Seidelmann, 1982; Kaplan, 1981). 

In this equatorial frame, some definition of the origin of right ascension must be made. We will 
not discuss that in this report, since one definition is at most a rotation from some other definition, 
and can be applied at any time. The important point is that consistent definitions must be used 
throughout the model development. The need for this consistency will, in all probability, eventually 
lead to our defining the origin of right ascension by means of the JPL planetary ephemerides, followed 
by our using interferometric observations of both natural radio sources and spacecraft at planetary 
encounters as a means of connecting the planetary and the radio reference frames (Folkner et al. } 
1994; Dewey, 1991; Newhall et aL f 1986). 

Also, unless otherwise stated, we will mean by “terrestrial reference frame” some reference frame 
tied to the mean surface features of the Earth. Currently, we are using a right-handed version of 
the Conventional International Origin (CIO) reference system with the pole defined by the 1903.0 
pole. In practice, this is accomplished by defining the position of one of the interferometric observing 
stations (generally Deep Space Station (DSS) 15 at the Goldstone deep space tracking complex), and 
then by measuring the positions of the other stations under a constraint. This constraint is that the 
determinations of Earth orientation agree on the average with the International Earth Rotation Service 
(IERS) (1993) |and its predecessor, Bureau International de l’Heure (BIH) (1983)] measurements of 
the Earth’s orientation over some substantial time interval (« years). This procedure, or its functional 
equivalent, is necessary since the interferometer is sensitive only to the baseline vector. The VLBI 
technique does not have any preferred origin relative to the structure of the Earth. The rotation of 
the Earth does, however, provide a preferred direction in space which can be associated indirectly 
with the surface features of the Earth. 

In contrast, geodetic techniques which involve the use of artificial satellites, or the Moon, are 
sensitive to the center of mass of the Earth as well as the spin axis. Thus, those techniques require 
only a definition of the origin of longitude. Laser ranging to the retroreflectors on the Moon allows a 
realizable practical definition of a terrestrial frame, accurately positioned relative to a celestial frame 
which is tied to the planetary ephemerides (Folkner et a/,, 1994). The required collocation of the laser 
and VLBI stations is being provided by Global Positioning Satellite (GPS) measurements of baselines 
between VLBI and laser sites starting in the late 1980s ( e.g ., Ray et ai., 1991). Careful definitions 
and experiments of this sort are required to realize a coordinate system of centimeter accuracy. 

The relativistic delay formulation presented in this report is essentially identical to that in the 
original 1983 version (Fanselow, 1983). Special relativistic terms in the model delay have not been 
changed, but a small revision was made in the default gravitational correction starting with Rev. 
3 (Sovers and Fanselow, 1987), when TDT (Temps Dynamique Terrestre) spatial coordinates were 
introduced (see Sec. 2.1.3). Among the estimated parameters, only baseline lengths were affected by 
this change, in that all distances were increased by the same factor of % 2 parts in 10 s . 

Except for subcentimeter relativistic complications caused by the locally varying Earth potential 
(as discussed below), calculation of the VLBI model for the observed delay can be summarized as: 
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1. Specify the proper locations of the two stations as measured in an Earth-fixed frame at the time 
that the wave front intersects station #1. Let this time be the proper time t' x as measured by a 
clock in the Earth-fixed frame. 

2. Modify the station locations for Earth-fixed effects such as solid Earth tides, tectonic motion, 
and other local station motion. 

3. Transform these proper station locations to a geocentric celestial coordinate system with its 
origin at the center of the Earth, but moving with the Earth. This is a composite of 12 separate 
rotations, represented by a rotation matrix Q(t). 

4. Perform a Lorentz transformation of these proper station locations from the geocentric celestial 
frame to a frame at rest relative to the center of mass of the Solar System, and rotationally 
aligned with the celestial geocentric frame. 

5. In this Solar-System-barycentric (SSB) frame, compute the proper time delay for the passage of 
the specified wave front from station #1 to station #2. Correct for source structure. Add the 
effective change in proper delay caused by the differential gravitational retardation of the signal. 

6. Perform a Lorentz transformation of this SSB geometric delay back to the celestial geocentric 
frame moving with the Earth. This produces the adopted model for the geometric portion of the 
observed delay. 

7. To this geometric delay, add the contributions due to clock offsets, to tropospheric delays, and 
to the effects of the ionosphere on the signal (see Sections 3 through 5). 

As indicated in step 5, the initial calculation of delay is carried out in a frame at rest relative to 
the center of mass of the Solar System (SSB frame.) First, however, steps 1 through 4 are carried out 
in order to relate proper locations in the Earth-fixed frame to corresponding proper locations in the 
SSB frame. Step 4 in this process transforms station locations from the geocentric celestial frame to 
the SSB frame. This step incorporates special-relativistic effects to all orders of the velocity ratio v/c. 
In the presence of gravity, this transformation can be viewed as a special relativistic transformation 
between proper coordinates of two local frames (geocentric and SSB) in relative motion. For both 
frames, the underlying gravitational potential can be viewed approximately as the sum of locally 
constant potentials caused by all masses in the Solar System. The complications caused by small local 
variations in the Earth’s potential are discussed below. Initial proper delay is then computed (step 5) 
in the SSB frame on the basis of these SSB station locations and an a priori SSB source location. A 
small proper-delay correction is then applied to account for the differential gravitational retardation 
introduced along the two ray paths through the Solar System, including retardation by the Earth’s 
gravity. A final Lorentz transformation including all orders of v/c then transforms the corrected SSB 
proper delay to a model for the observed delay. 

Since the Earth’s potential varies slightly across the Earth (At/j^/c 2 « 4 x 10“ 10 from center 
to surface), the specification of proper distance is not as straightforward with respect to the Earth’s 
potential as it is with respect to the essentially constant potentials of distant masses. To overcome this 
difficulty, output station locations are specified in terms of the “TDT spatial coordinates” (Shahid- 
Saless tt al.y 1991) used in Earth-orbiter models. Baselines modeled on the basis of this convention 
deviate slightly in length (< 2 mm) from the proper values. A proper length that corresponds to a 
modeled baseline can be obtained through appropriate integration of the local metric (Shahid-Saless 
tt a/., 1991). In practice, such a conversion is not necessary since comparison of baseline measurements 
obtained by different groups would be carried out in terms of TDT spatial coordinates. 

The current model has been compared (Thomas, 1991; TYeuhaft, 1991) with the “1-picosecond” 
relativistic model for VLBI delays developed by Shahid-Saless tt aL (1991). When reduced to the same 
form, the model presented here is identical to that model at the picosecond level, term by term, with 
one exception. Treuhaft and Thomas (1991) show that a correction is needed to the Shahid-Saless 
tt al. SSB system modeling of the atmospheric delay. This correction changes the Shahid-Saless tt aL 
result by as much as 10 picoseconds. The remainder of this section provides the details for the first 
six steps of the general outline above. 
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2.1 TIME INTERVAL FOR THE PASSAGE OF A WAVE FRONT 
BETWEEN TWO STATIONS 

The fundamental part of the geometric model is the calculation (step #5 above) of the time 
interval for the passage of a wave front from station to station $2. This calculation is actually 
performed in a coordinate frame at rest relative to the center of mass of the Solar System. This part 
of the model is presented first to provide a context for the subsequent sections, all of which are heavily 
involved with the details of time definitions and coordinate transformations. We will use the same 
subscript and superscript notation which is used in Section 2.7 to refer to the station locations as seen 
by an observer at rest relative to the center of mass of the Solar System. 

First, we calculate the proper time delay that would be observed if the wave front were planar. 
Next, we generalise this calculation to a curved wave front, and finally, we take into account the 
incremental effect which results from the fact that we must consider wave fronts that propagate 
through the various gravitational potential wells in the Solar System. 

2.1.1 PLANE WAVE FRONT 



Figure 1* Geometry for calculating the transit time of a plane wave front 


Consider the case of a plane wave moving in the direction, k, with station 2 having a mean 
velocity, as shown in Figure 1. As mentioned above, distance and time are to be represented as 
proper coordinates in the SSB frame. The speed of light c is set equal to 1 in the following formulation. 
The proper time delay is the time it takes the wave front to move the distance l at speed c. This 
distance is the sum of the two solid lines perpendicular to the wave front in Figure 1: 

*2 -ti = k ■ [r 2 (*i) - rx(ti)] + k ■ fi 2 \t ‘ 2 - ti] (2.1) 

where the superscript * serves to emphasize that station #2 has moved since tj 
following expression for the geometric delay: 

. £ NUi) -ri(ti)) 


. This leads to the 
( 2 . 2 ) 
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The baseline vector, r 2 (tx) - nfa), is computed on the basis of proper station locations calculated 
according to Eq. (2.183) in Sec. 2.7, 


2.1.2 CURVED WAVE FRONT 

In the case of a signal generated by a radio source within the Solar System it is necessary to 
include the effect of the curvature of the wave front. As depicted in Figure 2, let a source irradiate two 
who,, position* piv.n by „(fl r.L.iv. to th« Earth'. c,nt„. Th« pos.t.on o 
the Earth’s center, R e (ti), as a function of signal reception tune, t u at statwn #1 is mewure e 
to the position of the emitter at the time, of emission of the signal received at time ‘^Whilethis 
calculation is actually done in the Solar System barycentric coordinate system, the development that 
follows is by no means restricted in applicability to that frame. 



Figure 2. Geometry for calculating the transit time of a curved wave front 


Suppose that a wave front emitted by the source at time t e reaches station #1 at time tj and 
arrives at station #2 at time t' 2 . The geometric delay in this frame will be given by: 

T = t;- tl = |R 2 (t*)|-|Ri(ti)| ( 2 - 3 ) 

where all distances are again measured in units of light travel time. If we approximate the velocity of 
station #2 by 


£2 = 


R 2 (G) -R 2 (ti) 


and use the relation 


t\ - fi 

R,-(ti) = Re(ti) +*(<1) 


(2.4) 
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we obtain: 


r = |R c (*i) + r 2 (*i) + 02 t \ “ + ri x) 1 

= -Rc(^l) [ | fee + e 2 \ “ | fee + *l] ] 


( 2 . 6 ) 


where 


and 


£ 2 (^ 1 ) + ^ 2 r 
R c (t 1 ) 



(2.7) 

( 2 . 8 ) 


For ei and e 2 < 10 4 , we need to keep only terms of order e 3 in a sixteen-place machine in order to 
expand the expression for r in equation (2.6). This gives us: 


f + fi cA c (r) 

[1-ftc^a] 2 [l — R c ^ 2 ] ' ‘ ’ 

where to order c 3 

A«(r) = [el - el] - [(R c £ 2 ) 2 + (R c -«i) 2 + (R c e 2 ) 3 - (Rc « 2 )^ - (R c «i) 3 + (Rc«i)e?] (2.10) 

The first term in (2.9) is just the plane wave approximation , i.e., as R c — ♦ oo, R c — * k, with the 
second term in brackets in (2.10) approaching zero as r 2 /iZ c . Given that the ratio of the first term 
to the second term is « r/R c , wave front curvature is not calculable in a sixteen-place machine for 
R > 10 16 x r. For Earth-fixed baselines that are as long as an Earth diameter, requiring that the 
effects of curvature be less than 0.01 cm implies that the above formulation (2,10) must be used for 
R < 1.4 x 10 15 km, or approximately 150 light years. 

The procedure for the solution of (2.9) is iterative for e < 10“ 4 , using the following: 


Tn 


-RcA c (r n _i) 

To -f- 

2 ( 1 — 


where 


*0 — Tpi ane wave 


For € > 10~ 4 , directly iterate on the equation (2.6) itself, using the procedure: 


( 2 . 11 ) 


( 2 . 12 ) 


T n — #c|Rc + *2( r n-l)| “ i2c|R C + *l| 


(2.13) 


where again tq is the plane wave approximation. 


2.1.3 GRAVITATIONAL DELAY 


Because a light signal propagating in a gravitational potential is retarded relative to its motion 
in field-free space, the computed value for the differential time of arrival of the signals at rj(ti) and 
r 2 (i 2 ) must b e corrected for gravitational effects. Gravitational potential effects and curved wave 
front effects are calculated independently of each other since the former are a small perturbation 
(» 8.5 microradians or < 1. 75), even for Sun-grazing rays. 

For the geometry illustrated in Figure 3, the required correction to coordinate time delay is given 
by Moyer (1971) as: 


(1 + Irrs'lUv 

In 

r, -f r 2 (t 2 ) + r $ 2 

— In 

r. + ri(tj) + r,i 

C 3 

111 

r„ -r T 2 (t' 2 ) - t, 2 

r, + - t, : 


(2.14) 
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where r # , is defined as: 


r*i = |r i(U) - r,(t e )| 

Here 7 rrhr is the 7 factor in the parametrized post-Newtonian gravitational theory (e.g 
1973). In the Brans-Dicke theory, 

_ 1 + w 

1 ” N ~ 2+Z 

where w is the coupling constant of the scalar field. For general relativity, = 1, i.e., w — ► 00 . 

However, we allow to be an estimated parameter so that by setting 7 j*fn — — 1, we also have 

the option of “turning off” the effects of general relativity on the estimate of the delay. This proves 
useful for software development. The gravitational constant, /x p , is 

= Cm? (2.17) 

where G is the universal gravitational constant, and rrip is the mass of the pth gravitating body. 


(2.15) 
. Misner ct al. t 

(2.16) 



Figure S. A schematic representation of the geodesic connecting two points in the 
presence of a gravitational mass 


Depending on the particular source-receiver geometry in a VLBI experiment, a number of ap- 
proximations are possible for the correction Aq p of (2.14). Dropping the time arguments in (2.14), 
we have: _ _ _ _ 


A g p — 


( 1 )Mp 


r, + r 2 + r s7 


r, + ri - r,i 

r, +r, + r„_ 


r. + r 2 - r,2 J 


(2.18) 


This formulation is fine for r s R£ r t « r 9X , but can be put in a computationally better form for the 
case of distant sources with closely spaced VLBI receivers, i.e., |r 2 — rj|/ri — * 0, r\/r, — 0. For these 
sources, expand A g p in terms of r, jr 9 , r Jt /r,, and make use of the relationsliip 


= r; 


_ 21 1 / 2 ~ 

- 2r„ - r, + rf] « r, - r, ■ r. 


(2.19) 
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This leads to 


^Gp 


(I ~h 'Tffy ) ^ 


r i iTrfi 

r 3 -f r 2 r. 


for r*/r, — ► 0. 

If we further require that |r 2 — rxl/r; — • 0, and make use of 


12 = 11 + Ar 


then: 

r 2 +r 2 r, = rj ^1 + 2 ?i ■ Ar/rx + (Ar/r!) 2 j / +ri f, + Ar f, 
ss r, (l + fi • Ar /r,) + ri r, 4- Ar • r. 

In the limit of A t/t\ —* 0: 

r 2 (l + r 2 • T t ) — * ri(l +- ri r,) + Ar * (?i + r ( ) 

Substituting into (2.20) and expanding the logarithm, we obtain: 

(A )m p (r 2 — r x ) (r x -+r,) 


Ag p = — - 


r i(l +?i r,) 


( 2 . 20 ) 


( 2 . 21 ) 


( 2 . 22 ) 


(2.23) 


(2.24) 


Using whichever of these three formulations (2.18, 2.20 or 2.24) is computationally appropriate, 
the model calculates a correction Ac P for each of the major bodies in the Solar System (Sun, planets, 
Earth, and Moon). 

Before the correction Ac?p can be applied to a proper delay computed according to Eq. (2.2), it 
must be converted from a coordinate-delay correction to a proper-delay correction appropriate to a 
near-Earth frame. For such proper delays, the gravitational correction is given to good approximation 
by 

A'g p = Ag p -(Hw)^ (2.25) 

where r is the proper delay given by Eq. (2.2), and where U is the negative of the gravitational 
potential of the given mass divided by c 2 , as observed in the vicinity of the Earth (U is a positive 
quantity). The Ur term is a consequence of the relationship of coordinate time to proper time, and 
the term is a consequence of the relationship of coordinate distance to proper distance. 

The total gravitational correction used is: 


N 


A G — XI ^ 


Gp 


P=1 


(2.26) 


where the summation to N is over the major bodies in the Solar System. For the Earth, the 
(1 +■ IrrsW 7 term in Eq. (2.25) is omitted if one wishes to conform with the “TDT spatial coordi- 
nates” used to reduce Earth-orbiter data. The scale factor (l + 7 Ff **)U is approximately 1.97 x 10“ 8 
for the Sun. A number of other conventions are possible. One of these, which does not omit the 
(1 +• ') rrff )Ur term for the Earth, but evaluates it at the Earth’s surface, yields an additional scale 
factor of 0.14 x 10“ 8 . In either case, the model delay is decreased. Consequently, all inferred “mea- 
sured” lengths increase by the same fraction relative to previous lengths (c.g. by 19.7 parts per billion 
or 21.1 ppb). 

Some care must be taken in defining the positions given by r,, and ri(tx)- We have chosen 

as the origin the position of the gravitational mass at the time of closest approach of the received 
signal to that object. The position, r,, of the source relative to this origin is the position of that 
source at the time, £«, of the emission of the received signal. Likewise, the position, r,(£ t ), of the ith 
receiver is its position in this coordinate system at the time of reception of the signal. Even with 
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this care in the definition of the relative positions, we are making an approximation, and implicitly 
assuming that such an approximation is no worse than the approximations used by Moyer (1971) to 
obtain (2.14). 

Some considerations follow, regarding the use of appropriate times to obtain the positions of the 
emitter, the gravitational object, and the receivers. For a grazing ray emitted by a source at infinity, 
using the position of the gravitating body G at the time of reception of the signal at station #1 rather 
than at the time of closest approach of the signal to G can cause a 15-cm error on baselines with 
a length of one Earth radius as shown by the following calculation. From Figure 4, the calculated 
distance of closest approach, R> changes during the light transit time, tught tranaiti of a signal from a 
gravitational object at a distance Reg by: 

A R « Reg © * tught tran§it — © ' R%g/ c (2.27) 


EMITTER 



Figure 4. A schematic representation of the motion of a gravitating object during the 
transit time of a signal from the point of closest approach to reception by 
an antenna 


Since the deflection is: 


A© a 2 


(i-’-'WkK. 



(2.28) 
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S{ A©) : 


[— 1 

= Ae[®*M 

= Ae\ REG ] 

| 50 

l R J 

Lc^ceJ 

L c0 J 

1 dt 


(2.29) 


We consider the two bodies of largest mass in the Solar System: the Sun and Jupiter. For grazing 
rays, their respective deflections A0 are 8480 and 73 nanoradians. The barycentric angular velocities 
are estimated to be 0.06 and 17 nrad/s for the Sun and Jupiter. Note that Eq. (2.27) does not 
apply to the Sun. The Sun’s motion in the barycentric frame has a period of 11 years with a radius 
on the order of the Sun’s radius. Using approximate radii and distances from Earth to estimate Reg 
and 0, Eq. (2.29) gives 25 nrad for Jupiter; the corresponding value for the Sun is 0.07 nrad. For 
a baseline whose length equals the radius of the Earth, £(A0).R£ is thus approximately 0.05 and 
15 cm for the Sun and Jupiter, respectively. The effect is much smaller for the Sun in spite of its 
much larger mass, due to its extremely slow motion in the barycentric frame. 

In view of the rapid decrease of gravitational deflection with increasing distance of closest ap- 
proach, it is extremely unlikely that a routine VLBI observation would involve rays passing close 
enough to a gravitating body for this correction to be of importance. Exceptions are experiments 
specifically designed to measure planetary gravitational bending (Treuhaft and Lowe, 1991). In order 
to guard against such an unlikely situation in routine work, and to provide analysis capability for spe- 
cial experiments, the MODEST code always performs the transit-time correction for all planets. To 
obtain the positions of the gravitational objects, we employ an iterative procedure, using the positions 
and velocities of the objects at signal reception time. If R(t r ) is the position of the gravitational object 
at signal reception time, t P , then that object’s position, R(t a ), at the time, t OJ of closest approach of 
the ray path to the object was: 

R(t 0 ) = R(t r ) - Y\t r - t„l (2.30) 


t r t a — 


|Ro 


(2.31) 


We do this correction iteratively, using the velocity, V(t r ), as an approximation of the mean velocity, 
V. Because v/c « 10~ 4 , an iterative solution, 


^Ln(^o) — R(£r) 


V(M 


|Rn-x(t.)| 


(2.32) 


rapidly converges to the required accuracy. 
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2.2 TIME INFORMATION 

Before continuing the description of the geometric model, a few words must be said about 
time-tag information and the time units which will appear as arguments below. A general reference 
for time definitions is the Explanatory Supplement, 1992. The epoch timing inform; at.on n he c data 
is taken from the UTC (Universal Coordinated Time) time tags in the data stream at station #1, 
UTC, TiTtim, is converted to Dynnn.it Time [TDT) and is aUo need ns nn „ 8 „ m en. .0 

obtain an a priori estimate of Earth orientation. The conversion consists of the following components. 


TDT = (TDT- TAI ) -f (TAJ - UTCjers ) + ( UTCiers ~ UTCstd) 
+ ( UTCstd - UTC X ) + tfTC, 


(2.33) 


The four offsets in (2.33) thus serve to convert the station 1 time tags to TDT. In turn, their meaning 
is the following: 

1 TDT -TAI is 32.184 seconds by definition; TAI (Temps Atomique International) is atomic time, 
o TAI - UTC if ns is the offset between atomic and coordinated tune. The International Earth 
Rotation Service (IERS), its predecessor, Bureau International de l’Heure (BIH), and Bureau 
International des Poids et Mesures (BIPM) are the coordinating bodies responsible for upkeep 
and publication of standard time and Earth rotation quantities. TAI - UTCiers « a published 
integer second offset (leap seconds) for any epoch after 1 January, W72. that time. i lb 
a more complicated function, which need not be discussed here since MODEST is only intended 

to model observations starting in the mid-1970s. , , , , 

3. UTCiers ~ UTCstd is the offset in UTC between the coordinated time scales maintained by 
the IERS (BIPM) and secondary standards maintained by numerous national organizations^ For 
VLBI stations in the U.S. this secondary standard is that of the National Institute of Standards 
and Technology (NIST) in Boulder, Colorado. These offsets can be obtained from BIPM Circular 

T (e.g. BIPM, 1990). . , . . . 

4. UTC - UTC\ is the (unknown) offset between UTC kept by station 1 and the secondary na- 
tional standard. This normally amounts to several /is, but may not be precisely known for 
each experiment. It is a source of modeling error: an error St in epoch time causes an error of 
ss Bu E At = 7.3 x 10“ 6 cm per km baseline per fis of clock error, where u E is the rotation rate 
of the Earth (Section 7). For the extreme case of a 10,000 km baseline, however, this amounts to 
only 0.07 cm per /is clock offset. 

A priori UTl-UTC and pole positions are normally obtained by interpolation of the IERS 
Bulletin A smoothed values. However, any other source of UTl-UTC and pole position could be 
used provided it is a function of UTC, and is expressed in a left-handed coordinate system (see 
Section 2.6.1). Part of the documentation for any particular set of results should clearly state what 
were the values of UTl-UTC and pole position used in the data reduction process. 

For the Earth model based on the new IAU conventions, the following definitions are employed 

throughout (Kaplan, 1981): 

1. Julian date at epoch J2000 = 2451545.0. _ . , 

2. All time arguments denoted by T below are measured in Julian centuries of 36525 days of the 
appropriate time relative to the epoch J2000, i.e., T = ( JD - 2451545.0)/36525. 

3 For the time arguments used to obtain precession, nutation, or to reference the ephemens 
Barycentric Dynamic Time ( TDB , Temps Dynamique Barycentrique) is used. This is related 
to Terrestrial Dynamic Time (TDT, Temps Dynamique Terrestre) by the following: 

TDB = TDT + 0. ‘001658 sin(g + 0.0167 sin(g)) ( 2 - 34 ) 


where 


(357,°528 + 35999. °050 TDT] x 2r 
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2.3 STATION LOCATIONS 

Coordinates of the observing stations are expressed in the Conventional International Origin 
(CIO) 1903.0 reference system, with the reference point for each antenna defined as in Sec. 2.8. 
The pre-1984 model considered the three coordinates of station i: r #p ., A*, z, (radius off spin axis, 
longitude, and height above the equator, respectively) to be time- invariant. At the present accuracy 
level of space geodesy, it is imperative to account for tectonic motion. This is most simply done by 
estimating a new set of coordinates in the least-squares process for each VLBI session. Post-processing 
software then makes linear fits to these results to infer the time rate of change of the station location. 
Care must be taken that the correlations of coordinates estimated at different epochs are accounted 
for properly. The advantage of this approach is that the contribution of each session to the overall 
slope may be independently evaluated, since it is clearly isolated. An alternative second approach is 
to model tectonic motion directly, and to introduce time rates of change of the station coordinates as 
parameters in MODEST. The model is linear, with the cylindrical coordinates at time t expressed as 

r >Pi = r *p, + r.p, (* - to) (2-36) 

= A? + A f (t-t 0 ) (2-37) 

*«*?+*(* -*o) (2-38) 

Here t 0 is a reference epoch, at which the station coordinates are (r° p ., A°, r°). If modeling is done 
in Cartesian coordinates, the analogous expressions are 

*i = *? + ±i(t- to) (2.39) 

i u = y° + y «i(f — t 0 ) (2.40) 

Zi = z° + i<(t “ to) (2-41) 

with (xP,t/P,zP) being the station coordinates at the reference epoch to * 

2.3.1 MODELS OF TECTONIC PLATE MOTION 

As alternatives to estimating linear time dependence of the station coordinates, several 
standard models of tectonic plate motion are optionally available in MODEST. They all describe the 
motion as a rotation of a given plate about its rotation pole on the surface of a spherical Earth. Time 
dependence of the Cartesian station coordinates of station i which resides on plate j is expressed as 


Xi = + (wjz,° - wjy°)(t - to) 

(2.42) 

y. = y? + (^i 1 ? - w i*. 9 )(* - *o) 

(2.43) 

*-^ + («^l£-«j*?)(t-t 0 ) 

(2.44) 


where a;£ e are the angular velocities. 

Although these models are based on paleomagnetic data spanning millions of years, they have 
been found to provide a good quantitative characterization of present-day plate motions. The model 
due to Minster and Jordan (1978) was the first to be used in VLBI analyses. It is denoted AMO-2 in 
the original paper, and is also described in an addition to the MERIT standards document (Melbourne 
et a/., 1983, 1985). Table I gives a list of the rotation rates for the 11 plates in the AMO-2 model. 
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Table I 

Plate Rotation Velocities: Minster- Jordan AMO-2 Model^ 


Plate 


Wy 

w x 

AFRC 

0.988 

-3.360 

4.192 

ANTA 

-0.923 

-1.657 

3.765 

ARAB 

4.867 

-2.922 

6.520 

CARB 

-0.486 

-0.988 

1.881 

COCO 

-11.122 

-23.238 

12.663 

EURA 

-0.536 

-2.769 

3.422 

INDI 

8.443 

4.365 

7.528 

NAZC 

-1.586 

-9.299 

11.006 

NOAM 

0.576 

-3.984 

-0.249 

PCFC 

-2.143 

5.439 

-11.438 

SOAM 

-0.978 

-1.863 

-1.508 


f units are nr ad /year 


Note that the velocities are expressed in nanoradians per year rather than the microdegrees per year 
used in the original paper. 

More recent models, denoted NUVEL-1 and NNR-NUVELl, are due to DeMets ct aL (1990) 
and Argus and Gordon (1991), respectively. In NUVEL-1, the Pacific plate is stationary, while 
NNR-NUVELl is based on the imposition of a no-net-rotation (NNR) condition. With some notable 
exceptions, the NUVEL models give rates that are very close to those of the AMO-2 model. The 
AMO-2 INDI plate has been split into AUST and INDI, and there are two additional plates: JDEF (Juan 
de Fuca) and PHIL (Philippine). The NUVEL-1 rotation rates are given in Tables II and III. 


Table II 

Plate Rotation Velocities: NUVEL-1 Model^ 


Plate 


Wy 

w x 

AFRC 

2.511 

-8.303 

14.529 

ANTA 

0.721 

-6.841 

14.302 

ARAB 

8.570 

-5.607 

17.496 

AUST 

9.777 

0.297 

16.997 

CARB 

1.393 

-8.602 

12.080 

COCO 

-9.323 

-27.657 

21.853 

EURA 

0.553 

-7.567 

13.724 

INDI 1 

8.555 

-5.020 

17.528 

JDEF 

6.81 

3.32 

5.31 

NAZC 

-0.023 

-14.032 

20.476 

NOAM 

1.849 

-8.826 

10.267 

PCFC 

0.000 

0.000 

0.000 

PHIL 

11.9 

12.8 

0.000 

SOAM 

0.494 

-6.646 

9.517 


t units are nrad/year 
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Table III 

Plate Rotation Velocities: NNR-NUVELl Model^ 


Plate 

w* 

U?y 

w z 

AFRC 

0.929 

-3.239 

4.098 

ANTA 

-0.862 

-1.777 

3.871 

ARAB 

6.987 

-0.543 

7.067 

AUST 

8.194 

5.362 

6.566 

GARB 

-0.190 

-3.538 

1.649 

COCO 

-10.907 

-22.592 

11.420 

EURA 

-1.030 

-2.503 

3.293 

INDI 

6.973 

0.045 

7.097 

JDEF 

5.227 

8.386 

-5.124 

NAZC 

-1.607 

-8.968 

10.046 

NOAM 

0.265 

-3.761 

-0.164 

PCFC 

-1.583 

5.065 

-10.430 

PHIL 

10.320 

-7.700 

-10.430 

SOAM 

-1.089 

-1.581 

-0.913 


f units are nrad/year 


A recent revision of the paleomagnetic time scale has led to a rescaling of the NUVEL rates. These 
“NUVEL-lA* and “NNR-NUVELlA” model rates are equal to the NUVEL-1 and NNR-NUVELl 
rates of Tables II and III, respectively, multiplied by a factor of 0.9562 (DeMets et al., 1994). Partial 
derivatives with respect to the plate velocities are given in Sec. 2.9.2. 
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2.4 TIDAL EFFECTS 


As an initial step in calculating the geometric delay, we need to consider the effects of crustal 
motions on station locations. Among these deformations are solid Earth tides, tectonic motions, and 
alterations of the Earth’s surface due to local geological, hydrological, and atmospheric processes. 
If only the crustal movement due to solid Earth tides and tectonic motion is modeled, then the 
remaining effects will manifest themselves as temporal changes of the Earth-fixed baseline. It is 
therefore important to model all crustal motions as completely as possible. The current level of 
mismodeling of these motions is probably one of the two biggest sources of systematic error (along 
with the troposphere) in analyses of VLBI data. 

In the standard terrestrial coordinate system, tidal effects modify the station location To by an 
amount 

A — A| 0 j 4" A po / &ocn &atm (2.45) 

where the four terms are due to solid Earth tides, pole tide, ocean loading, and atmosphere loading, 
respectively. Other Earth-fixed effects would be incorporated by augmenting the definition of A. 
All four tidal elTects are most easily calculated in some variant of the VEN (vertical, East, North) 
local geocentric coordinate system. To transform them to the Earth-fixed coordinate frame, the 
transformation VW t given in the next section, is applied. 

2.4.1 SOLID EARTH TIDES 

Calculating the alteration of the positions of the stations caused by solid Earth tides is rather 
complicated due to the solid tides’ coupling with the ocean tides, and the effects of local geology. 
We have chosen to gloss over these complications initially, and to incorporate the simple multipole 
response model described by Williams (1970), who used Melchior (1966) as a reference. Let R p be 
the position of a perturbing source in the terrestrial reference system, and ro the station position in 
the same coordinate system. To allow for a phase shift [rp) of the tidal effects from its nominal value 
of 0, the phase-shifted station vector r, is calculated from ro by applying a matrix X, describing a 
right-handed rotation through an angle t p about the Z axis of date, r, = Xro. This lag matrix, X, is: 

( cos rp sin rp 0 \ 

-sin rp cos rp 0 J (2.46) 

0 0 ij 

By a positive value of rp we mean that the peak response on an Earth meridian occurs at a time 
6t — rp/w E after that meridian containing r 0 crosses the tide-producing object, where u E is the 
angular rotation rate of the Earth. In the vertical component, the peak response occurs when the 
meridian containing r, also includes R p . 

The tidal potential at r, due to the perturbing source at R p is expressed as 

UU<ial = ^[(^) 2P2(C ° S ^ + (^) 3p 3( cos «)] 

= U 2 4- U 3 (2-47) 

where only the quadrupole and octupole terms have been retained. Here, G is the gravitational 
constant, m p is the mass of the perturbing source, P 2 and P 3 are Legendre polynomials, and 0 is the 
angle between r, and R p . While the quadrupole displacements are on the order of 50 cm, the mass 
and distance ratios of the Earth, Moon, and Sun limit the octupole terms to a few mm. The octupole 
terms are optionally included in the MODEST code, but partials with respect to the Love numbers 
are available only for the quadrupole terms. An estimate of the retardation correction (employing 
the position of the tide-producing mass at a time earlier than that of the observation by an amount 
equal to the light-travel time) shows that this correction is well below 1 mm, and can therefore be 
neglected. 
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6 is 


In a local geocentric VEN coordinate system on a spherical Earth, the tidal displacement vector 



(2.48) 


where the g £ *(t — 2,3) are the quadrupole and octupole displacements. The components of 6 are 
obtained from the tidal potential as 



(2.49) 

(2.50) 

(2.51) 


where h{(i — 2,3) are the vertical (quadrupole and octupole) Love numbers, Z*(i = 2,3) the corre- 
sponding horizontal Love numbers, and A, and are the station longitude and latitude, and g the 
acceleration due to gravity, 

g = Gm E /r 2 (2.52) 

Using the relation between terrestrial and celestial coordinates, 


cos & = sin 4> t sin 6 P + cos <f> s cos 6 P cos(A, - a p ) 


(2.53) 


with q p , S p the right ascension and declination of the perturbing body, and cxq the right ascension 
of Greenwich, some algebra produces the following expressions for the quadrupole and octupole com- 
ponents of 6 in terms of the coordinates of the station {x 9i y 9 ,z 4 ) and the tide-producing bodies 
[X p} Y PJ Z p ): 


4 ” - E 


4 ” -T.— 


R b 


4 ” 


-E 


3m p 


*p 5 


(r. R P ) 


r (r, Rp) 2 rjR 2 

R£ L 2 6 J 

( r * Rp)( s «ip ~ y.X p )/y/x 2 + yj 


x^JTyf z p - 


- /*' 7 (x.X P + y,Y P ) 

v x . + y 2 . 


(2.54) 

(2.55) 

(2.56) 


= E f£r( r ‘ ^ |5(r. • Rp) 2 - 3 r 2 R 2 p 

p 


92 ] = E f 5 (r. • Rp ) 2 - t]rA (x,Yp - y.XpJ/v^FT^F 

P p l 


„ (3 > - 
03 - / . 


3M P r 2 

2/?J 


5(r, Rp) 2 - t 2 R 2 


V^F+TF ^p - 


V^F+1F 


(x.X p + y.F p ) 


(2.57) 

(2.58) 

(2.59) 


where /i p is the ratio of the mass of the disturbing object, p, to the mass of the Earth, and 

Rp = [X P> Y P1 Z p ] (2.60) 

is the vector from the center of the Earth to that body. The summations are over tide-producing 
bodies, of which we include only the Sun and the Moon. 

The above formulation implicitly assumes that the Love numbers h t and l x are independent of 
the frequency of the tide-generating potential. Proper treatment entails a harmonic expansion of 
Eqs. (2.54)-(2.59) and using a different set of h t , /, for each frequency component. Presently, only 
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the first six largest nearly diurnal components are allowed to have Love numbers that differ from the 
standard values given in Sec. 7. Each harmonic term is denoted by its historical name, if it exists, and 
the Doodson code (IERS, 1992) (c.g., K\ and Doodson number = 165555). The Doodson notation 
classifies the tidal components according to increasing speed. The correction to the solid tidal radial 
displacement for the kth harmonic term at station s is given by 

Sh^ k = <5h£#fc(\/5/247r) 3 sin <j> 9 cos 4> $ sin(A, + B k ) (2-61) 

where 6h * is the difference between the nominal quadrupole ( h 2 ) Love number (0.609) and the fre- 
quency dependent Love number (Wahr, 1979), H k is the amplitude of the kth harmonic term in the 
tide generating expansion from Cartwright and Edden (1973), 4> 0 is the geocentric latitude of the the 
station, A, is the East longitude of the station and 6 k is the kth harmonic tide argument. The Love 
numbers and tidal amplitudes are listed in Table IV. 


Table IV 

Frequency Dependent Solid Earth Tide Parameters 


Component (£) 

|9 

H k (mm) 

V>i (166554) 

0.937 

3 

(165565) 

0.514 

50 

Ki (165555) 

0.520 

369 

(165545) 

0.526 

-7 

P, (163555) 

0.581 

-122 

Oi (145555) 

0.603 

-262 


These optional corrections yield (Naudet, 1994) additional purely vertical station displacements 


(in mm) of: 

(166554) 6g x = 0.37 sin 2<f> $ sin(A, + a G + 0 (2.62) 

(165565) 6gi = -1.84 sin2<£, sin(A, + a G - H) (2.63) 

Ki (165555) 6g x = -12.68 sin2<£, sin(A, + a G ) (2.64) 

(165545) 6g x = 0.24 sin 2<f> $ sin (A, H- qg + G) (2.65) 

P x (163555) 5^ = 1.32 sin 2<t>, sin[A, + a G - 2(fl + F - D)] (2.66) 

O x (145555) 6gi = 0.62 sin 2tf>, sin[A, +a G - 2(fl 4- F)] (2.67) 


where A,, olq are the station latitude and longitude and Greenwich RA, respectively. The astro- 
nomical arguments F, D, Cl (mean anomaly of the Sun, mean argument of the latitude of the Moon, 
mean elongation of the Moon from the Sun, and the mean longitude of the ascending lunar node) are 
defined in Sec. 2.6.2. These displacements are then summed and used as the first order correction 
to each station’s vertical displacement. Horizontal corrections are presently ignored. Note that the 
largest correction, the K x term, is identical to that already recommended in 1985 by the MERIT 
standards (Melbourne tt al. y 1985). 

To convert the locally referenced displacement, 6 , which is expressed in the VEN system, to 
the Earth-fixed frame, two rotations must be performed. The first, W % rotates by an angle, 4> 9 
(station geodetic latitude), about the y axis to an equatorial system. The second, V, rotates about 
the resultant 2 axis by angle, —A, (station longitude), to bring the displacements into the standard 
geocentric coordinate system. The result is 

A 90 t = VWS ( 2 . 68 ) 
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where 


W = 


(2.69) 


and 



0 

1 

0 


— sin 
0 

COS 4 >* 



Actually, the product of these two matrices is coded: 


(2.70) 


( cos A, cos <f> a 
sin A, cos 4> $ 
sin <f> 9 


MODEST code uses geodetic latitudes 


4 >. 


— tan' 


->[- 

Lr . P , 


Z. 


.(1 - l//pj 


(2.71) 


(2.72) 


where /is the geoid flattening factor. The difference between geodetic and geocentric latitude can 
affect this model on the order of (tidal effect)/(flattening factor) « 0.1 cm. 


2.4.2 POLE TIDE 


One of the secondary tidal effects is the displacement of a station by the elastic response of 
the Earth’s crust to shifts in the spin axis orientation. The spin axis is known to describe a circle of 
« 20-m diameter at the north pole. Depending on where the spin axis pierces the crust at the instant 
of a VLBI measurement, the “pole tide* displacement will vary from time to time. This effect must 
be included if centimeter accuracy is desired. 

Yoder (1984) and Wahr (1985) derived an expression for the displacement of a point at geocentric 
latitude <f>,, longitude A, due to the pole tide: 


uiR, 


+ cos 2 cos A, + y sin A, ) 14> 
+ sin <j>, (— x sin A, + y cos A, ) / A] 


(2.73) 


Here ui E is the rotation rate of the Earth, R the radius of the (spherical) Earth, g the acceleration 
due to gravity at the Earth’s surface, and h and l the customary Love numbers. Displacements of the 
instantaneous spin axis from the current average spin axis along the x and y axes are given by x and 
y. Eq. (2.73)^shows how these map into station displacements along the unit vectors in the radial 
(r), latitude ( 4> ), and longitude (A) directions. With the standard values lj e = 7.292 x 10“ 5 rad/sec, 
R = 6378 km, and g = 980.665 cm/sec 2 , the factor u%R/g = 3.459 x 10“ 3 . Since the maximum values 
of * and y are on the order of 10 meters, and h « 0.6, l « 0.08, the maximum displacement due to 
the pole tide is 1 to 2 cm, depending on the location of the station ( <f> ,, A,). 

The locally referenced displacement 6 is transformed via the suitably modified transformation 
(2.71) to give the displacement A po i in the standard geocentric coordinate system. The pole tide 
effect has been coded as an optional part of the MODEST model. It is only applied if specifically 
requested, t.e., the default model contains no pole tide contributions to the station locations. 
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2.4.3 OCEAN LOADING 

This section is concerned with another of the secondary tidal effects, i.e., the elastic response 
of the Earth’s crust to ocean tides, which move the observing stations to the extent of a few cm. 
Such effects are commonly labeled “ocean loading." A model of ocean loading is incorporated in 
the MODEST code. It is general enough to accommodate a variety of externally derived parameters 
describing the tide phases and amplitudes at a number of frequencies. The present model entails 
deriving an expression for the locally referenced displacement 6 due to ocean loading. In a local 
Cartesian coordinate system (the computer code accepts inputs related to unit vectors in the Up, 
North, and West directions) at time t, 


6, =^^cos( W< t + V < -«/) (2-74) 

t=l 

The quantities u,< (frequency of tidal constituent t) and V) (astronomical argument of constituent i) 
depend only on the ephemeris information (positions of the Sun and Moon). The algorithm of Goad 
(IERS, 1989) is used to calculate these two quantities. On the other hand the amplitude and 
Greenwich phase lag Sf of each tidal component j are determined by the particular model assumed 
for the deformation of* the Earth. The local displacement vector is transformed via Eqs. (2.71) and 
(2.68) to the displacement A 0 cn in the standard geocentric frame. 

Input to MODEST provides for specification of up to 11 frequencies and astronomical arguments 
w, and V,, followed by tables of the local distortions and their phases, (j and 6j , calculated from the 
ocean tidal loading model of choice. The eleven components are denoted, in standard notation: K 2 , 
S 2 , A/ 2 , and N 2 (all with approximately 12-hour periods), Ki, Pi, 0\, Q\ (24 h), Mj (14 day), M m 
(monthly), and S, a (semiannual). 

Three choices of ocean loading models have been used with MODEST. They differ in the dis- 
placements calculated and components considered, as well as in the numerical values that they yield 
for the H s and */s. Scherneck’s results (1983, 1990, 1991) are the most complete in the sense of 
considering both vertical and horizontal displacements and all eleven tidal components. They have 
now been adopted for the IERS standards (1992), and compose the default ocean loading model. 
Goad’s model (1983) was adopted in the MERIT and early IERS standards (1989), but only considers 
vertical displacements. Pagiatakis’ (1982, 1990) model, baaed on Pagiatakis, Langley, and Vanicek 
(1982), considers only six tidal components (S 2 , A/ 2 , N 2 , Ki, Pi, and O j). 

An extension of the 1991 Scherneck model is also available in MODEST (Scherneck, 1993). The 
eleven tidal frequencies are modulated to an appreciable extent by multiples of N', the lunar nodal 
period (18.6 years). On the assumption that these additional terms yield ocean loading amplitudes 
which are in the same ratio to each main loading term as the companion tides are to the main tides, 
the additional station displacements can be written as 

r*i£l cos[(o;i + + v. + n ki N' - 6’)} (2-75) 

i=l k 

where the k summation extends over all integer multiples n k i of the lunar node N', and r kl is the 
ratio of the tidal amplitude of each companion k to the tidal amplitude of the parent t. Of 26 such 
components listed by Cartwright and Edden (1973), 20 are estimated to be significant in contributing 
to the largest ocean loading displacements at the 0.01 mm level. Table V shows the multiples n k> and 
amplitude ratios r^i for these 20 components. 
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Table V 

Lunar Node Companions to Ocean Tides 


1 

Component 

Companion 

r*»: Relative 
Amplitude 

K 2 (275555) 

-1 

-0.0128 


+ 1 

4-0.2980 


+2 

+0.0324 

S 2 (273555) 

-1 

4-0.0022 

M 2 (255555) 

-2 

4-0.0005 


-1 

-0.0373 

N 2 (245655) 

-i 

-0.0373 

K x (165555) 

-1 

-0,0198 


+1 

4-0.1356 


+2 

-0.0029 

Px (163555) 

-1 

-0.0112 

Ox (145555) 

-2 

-0.0058 


-1 

4-0.1885 

Qx (135655) 

—2 

4-0.0057 

j 

-1 

4-0.1884 

Mf (075555) 

+1 

4-0.4143 


+2 

4-0.0387 

M m (065455) 

-1 

-0.0657 


+1 

-0.0649 

S. a (057555) 

+1 

-0.0247 


In pushing the limits of Earth modeling to below 1 cm accuracy in the mid-1990s, ocean loading 
station displacements are one aspect of the models that are undergoing close scrutiny. Initial trials 
indicate that ocean loading amplitudes can be derived from VLBI experiments at an approximate 
accuracy level of 1-2 mm (Sovers, 1994). Partial derivatives of the VLBI observables with respect to 
eight diurnal and semidiurnal amplitudes are considered in Sec. 2.9.2. 1. When estimating parameters, 
however, great care must be used in order to avoid singularities due to the identity of components of 
station displacements (“confounding of parameters”). Since some components of ocean loading, solid 
Earth tides, and ocean tidally induced UTPM variations have the same frequencies, certain linear 
combinations of their station displacements are identical (see Sections 2.4.1 and 2.6. 1.1.2). 


2.4.4 ATMOSPHERE LOADING 

By analogy with the ocean tides that were considered in the previous section, a time-varying 
atmospheric pressure distribution can induce crustal deformation. A paper by Rabbel and Schuh 
( 1986 ) estimates the effects of atmospheric loading on VLBI baseline determinations, and concludes 
that they may amount to many millimeters of seasonal variation. In contrast to ocean tidal effects, 
analysis of the situation in the atmospheric case does not benefit from the presence of a well-understood 
periodic driving force. Otherwise, estimation of atmospheric loading via Green’s function techniques 
is analogous to methods used to calculate ocean loading effects. Rabbel and Schuh recommend a 
simplified form of the dependence of the vertical crust displacement on pressure distribution. It 
involves only the instantaneous pressure at the site in question, and an average pressure over a circular 
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region C of radius R = 2000 km surrounding the site. The expression for the vertical displacement 
(mm) is: 

A r = — 0.35p o - 0.55 p (2.76) 

where po is the local pressure anomaly (relative to the standard pressure of 1013.25 mbar), and p 
is the pressure anomaly within the 2000-km circular region mentioned above (both Quantities are in 
mbar). Note that the reference point for this displacement is the site location at standard (sea level) 
pressure. The locally referenced Ar is transformed to the standard geocentric coordinate system via 
the transformation (2.71). 

It was decided to incorporate this rudimentary model into MODEST as an optional part of the 
model, with an additional mechanism for characterizing p. The two-dimensional surface pressure 
distribution (relative to 1013.25 mbar) surrounding a site is described by 

p(x, y) = po 4- A\x 4- A 2 V + A$x 7 4- A 4 %y 4- A^xf 2 (2.77) 

where x and y are the local East and North distances of the point in question from the VLBI site. 
The pressure anomaly p may then be evaluated by the simple integration 

P = J J^dxdy p{x,y) / J J^dxdy (2.78) 

giving 

P = Po + (A 3 + A b )R 7 /4 (2.79) 

It remains the task of the data analyst to perform a quadratic fit to any available area weather data 
to determine the coefficients A1-5. Future advances in understanding the atmosphere-crust elastic 
interaction can probably be accommodated by adjusting the coefficients in Eq. (2.76). As an initial 
step along these lines, a station-dependent factor is introduced to scale the second coefficient in Eq. 
(2.76): 

Ar = — 0.35p o - 0.55(1 + f)p (2.80) 

This may account for differing geographical features surrounding different sites. In particular, / may 
depend on the fraction of ocean within the 2000 km radius. 

In summary, models have been presented that describe the four tidal effects on the stations 
(solid, pole, ocean, and atmosphere). Each of the locally referenced tidal displacement vectors is 
then transformed to the standard geocentric coordinate system via rotations like (2.68). After this 
transformation, the final station location is 

= r 0 4- A to t 4- A pot 4- Aocn 4" A a tm (2.81) 
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2.5 SOURCE STRUCTURE EFFECTS 


Numerous astrophysical studies during the past decade have shown that compact extragalactic 
radio sources exhibit structure on a milliarcsecond scale (e.p., Kellermann and Pauliny-Toth, 1981). 
Such studies are important for developing models of the origin of radio emission of these objects. 
Many radio source structures are found to be quite variable with frequency and time (Zensus and 
Pearson, 1987). If extragalactic sources are to serve as reference points in a stable reference frame, it 
is important to correct for the effects of their structures in astrometric VLBI observations. 

MODEST modeling allows optional corrections for the effects of source internal structures, based 
on work by Thomas (1980), Ulvestad (1988), and Chariot (1989, 1990a). A non-point- like distribution 
of the intensity of a source yields time dependent corrections to the group delay and delay rate 
observables, A r, and A f,, that may be written in terms of the intensity distribution J(e, u;,t) as 

At, = d<f> 9 /dw t At, = d<t> $ /dt (2.82) 

with 

<t>* = arctan(-Z,/Z c ) (2.83) 

and 

Z{ : } = / J dU /(., «, t) { “ } (2?rB «/A) (2.84) 

Here <fa 9 is the correction to the phase of the incoming signal, 8 is a vector from the adopted reference 
point to a point within the source intensity distribution in the plane of the sky, w and A are the 
observing frequency and wavelength, B the baseline vector, and the integration is over solid angles H. 
Source intensity distribution maps are most conveniently parametrized in terms of one of two models: 
superpositions of delta functions or Gaussians. At a given frequency, the corresponding intensity 
distributions are written as 

J(s) = ^2 s kH x - x k,y- yk) (2.85) 

k 

or 

c p 

•f( B ) = 2n -afcfc fc ex P [-K* “ *fc) cos & k + (y - yk) sin 6 k ] 2 / 2a* 2 

k 

- [-(x - X*) sin B k + (y - y fc ) cos fl*] 2 / 26 * 2 ] (2.86) 

where S k is the flux of component k , and Bk (with components x k) y k in the plane of the sky) is its 
position relative to the reference point. For Gaussian distributions, 8 * is the angle between the major 
axis of component k and the tt axis (to be defined below), and (0^,6*) are the full widths at half 
maximum of the (major, minor) axes of component k normalized by 2y/2log2. The quantities Z{ b j 
entering the structure phase 4>» [Eq. (2.83)] are 

8k / A ) ( 2 - 87 ) 

for delta functions, and 

*(;>=!> exp[-2* 2 (a 2 1/ 2 + fc 2 V* 2 )] { “ } (2*B s fc /A) (2.88) 

k 

for Gaussians. Here 

U k = u cos ^ 4- usinfl* 

V k = — it sin 8 k 4- v cos 6 k 
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( 2 . 89 ) 

( 2 . 90 ) 



with u, v being the projections of the baseline vector B on the plane of the sky in the E-W, N-S 
directions, respectively. 

MODEST accepts maps specified in terms of an arbitrary number of Gaussian or delta function 
components. At most, six parameters must be specified for each component: its polar coordinates 
and flux, and, for a Gaussian, its major and minor axes and the position angle of the major axis. The 
structural correction for phase is computed via Eqs. (2.83), (2.87), and (2.88). For the bandwidth 
synthesis (BWS) delay observable, the structure correction is the slope of a straight line fitted to 
the individual structure phases calculated for each frequency channel used during the observation. 
For example, for Mark III data there are typically 8 channels spanning «8.2 to 8.6 GHz at X band, 
and 6 channels spanning «2.2 to 2.3 GHz at S band. Delay rate structure corrections are calculated 
by differencing the structure phases at ±2 seconds (see Section 6). In the case of dual-band (S-X) 
experiments, a linear combination of the structure corrections calculated independently for each band 
is applied to the dual-band observables. 

The practical question to be resolved is whether such structural corrections based on maps yield 
significant and detectable corrections to the observables at the present levels of experimental and 
modeling uncertainty. Maps are available for only a few of the hundreds of sources currently observed 
by VLBI. Some of the extended sources show time variability on a scale of months; since the corrections 
At # and At, are quite sensitive to fine details of the structure, in such cases new maps may be required 
on short time scales. Depending on the relative orientation of the source and baseline, the delay 
correction can be as large as «1 ns, which is equivalent to tens of cm. Nevertheless, the prognosis 
appears to be good. Chariot (1990b) found that data from a multiple baseline geodynamics experiment 
are adequate to map source structures with high angular resolution. More recently Chariot (1993) 
has also shown that use of maps for the structure of the source 3C 273 improves the fit in analyses of 
geodetic experiments. 

Empirical evaluation of the effects of unknown source structure on VLBI measurements could be 
made via the time rates of change of the source right ascension a and declination 6. A linear model 
of the motion of source coordinates 

a = ao d(t — to) (2.91) 

6 = 6o + £(t * to) (2.92) 

is implemented in MODEST. Non-zero estimates of the rate parameters d and 6 could arise either 
from genuine proper motion or from motion of the effective source centroid sampled by VLBI mea- 
surements. Proper interpretation of such results is problematic, but non-zero rates can be used as a 
crude diagnostic for the presence of structure effects. 

MODEST code also provides the option of modeling source structure as a superposition of two 
6 functions centered at points Pi(xi,yi) and P 2 ( i 2 iV 2 ) respectively, as in Eq. (2.85) above. The 
parameters describing the two components are: 1) flux ratio K — S 2 /S 1 ) where Sk is the flux of the 

jfcth component, 2) component separation s = |s| = IP 1 P 2 I, and 3) position angle 6. The posit ion angle 

is 6 = 0° when P\P? is in the direction of increasing declination 6 , and $ = 90° when P 1 P 2 is in the 

direction of increasing right ascension d. From Chariot (1990b), the group delay has the following 
dependence on the structural parameters: 

2 *K (1 - K) R\l - cos(2t tR)} . Q3 % 

T w (1 + K) ' \K 2 + 2Kcos{2irR) + 1] 


where 


R = B s/A 


(2.94) 


For evaluating partial derivatives of r, the component separation s and baseline B are most conve- 
niently written in terms of their components in the celestial system, as 


Then R becomes 


s = ds sin# + 6 s cos 6 

(2.95) 

B = d u A + £ v A 

(2.96) 

R = s (u sin 0 + v cos#) / A 

(2.97) 
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2.6 TRANSFORMATION FROM TERRESTRIAL TO CELESTIAL 
COORDINATE SYSTEMS 

The Earth is approximately an oblate spheroid, spinning in the presence of two massive moving 
objects (the Sun and the Moon) which are positioned such that their time-varying gravitational effects 
not only produce tides on the Earth, but also subject it to torques. In addition, the Earth is covered 
by a complicated fluid layer, and also is not perfectly rigid internally. As a result, the orientation 
of the Earth is a very complicated function of time, which to first order can be represented as the 
composite of a time-varying rotation rate, a wobble, a nutation, and a precession. The exchange of 
angular momentum between the solid Earth and the fluids on its surface is not readily predictable, 
and thus must be continually determined experimentally. Nutation and precession are well modeled 
theoretically. However, at the accuracy with which VLBI can determine baseline vectors, even these 
models are not completely adequate. 

Currently, the rotational transformation, Q, of coordinate frames from the terrestrial frame to 
the celestial geocentric frame is composed of 6 separate rotations (actually 12, since the nutation, 
precession, and “perturbation” transformations, N } P } and fi, consist of 3 transformations each) 
applied to a vector in the terrestrial system: 

Q = 0 PNUXY (2.98) 

In order of appearance in (2.98), the transformations are: the perturbation rotation, precession, 
nutation, UTl, and the x and y components of polar motion. All are discussed in detail in the 
following four sections. With this definition of Q, if r t is a station location expressed in the terrestrial 
system, e.p., the result of (2.81), that location, r c , expressed in the celestial system is 


t c = Qr t (2.99) 

This particular formulation follows the historical path of astrometry, and is couched in that 
language. While esthetically unsatisfactory with modern measurement techniques, such a formulation 
is currently practical for intercomparison of techniques and for effecting a smooth inclusion of the 
interferometer data into the long historical record of astrometric data. Much more pleasing esthetically 
would be the separation of Q into two rotation matrices: 

Q = ( 2 . 100 ) 

where Q 2 are those rotations to which the Earth would be subjected if all external torques were 
removed (approximately UXY above), and where Qi are those rotations arising from external torques 
(approximately QPN above). Even then, the tidal response of the Earth prevents such a separation 
from being perfectly realized. Eventually, the entire problem of obtaining the matrix Q, and the tidal 
effects on station locations should be done numerically. Note that the six rotations operating on a 
vector yield its components in a new coordinate system, and, since we rotate the Earth rather than the 
celestial sphere, the matrices Q, P ) and N will be the transposes of those used to rotate the celestial 
system of J2000 to a celestial system of date. 

2.6.1 UTl AND POLAR MOTION 

The first transformation, Y, is a right-handed rotation about the x axis of the terrestrial frame 
by an angle 02- Currently, the terrestrial frame is the 1903.0 CIO frame, except that the positive y 
axis is at 90 degrees east (Moscow). The x axis is coincident with the 1903.0 meridian of Greenwich, 
and the z axis is the 1903.0 standard pole. 

f 1 

Y = 0 

Vo 


0 0 \ 

cos ©2 sin ©2 I (2.101) 

— sin ©2 cos ©2 J 
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where 02 is the y pole position published by IERS. 

The next rotation in sequence is the right-handed rotation (through an angle ©i about the y 
axis) obtained after the previous rotation has been applied: 

( cos ©i 0 — sin ©i \ 

01 0 ( 2 . 102 ) 
sin ©i 0 cos ©i J 

In this rotation, ©i is the IERS x pole position. Note that we have incorporated in the matrix 
definitions the transformation from the left-handed system used by IERS to the right-handed system 
we use. Note also that instead of IERS data used as a pole definition, we could instead use any other 
source of polar motion data provided it was represented in a left-handed system. The only effect 
would be a change in the definition of the terrestrial reference system. 

The application of U XY ” to a vector in the terrestrial system of coordinates expresses that vector 
as it would be observed in a coordinate frame whose z axis was along the Earth’s ephemeris pole. 
The third rotation, U , is about the resultant z axis obtained by applying U XY*. It is a rotation 
through the angle, -if, where H is the hour angle of the true equinox of date (t.e., the dihedral angle 
measured westward between the xz plane defined above and the meridian plane containing the true 
equinox of date). The equinox of date is the point defined on the celestial equator by the intersection 
of the mean ecliptic with that equator. It is that intersection where the mean ecliptic rises from below 
the equator to above it (ascending node). 


f cos H — sin H 0\ 

U = J sin H cos if 0 1 (2.103) 

V ° o i; 

This angle H is composed of two parts: 

if = h 7 + a £ (2.104) 

where h-, is the hour angle of the mean equinox of date, and as (equation of equinoxes) is the difference 
in hour angle of the true equinox of date and the mean equinox of date, a difference which is due 
to the nutation of the Earth. This set of definitions is cumbersome and couples the nutation and 
precession effects into Earth rotation measurements. However, in order to provide a direct estimate 
of conventional UT 1 (universal time) it is convenient to endure this historical approach, at least for 
the near future. 

UT 1 is defined to be such that the hour angle of the mean equinox of date is given by the following 
expression (Aoki tt al. t 1982; Kaplan, 1981): 


h , = UT 1 + 
+ 

6 h 41 m 50*. 54841 + 8640184*. 8 12866 T u 
0*. 093 104 Tl - 6* .2 x 10" 6 T„ 

(2.105) 

where the dimensionless quantity 

T 

(Julian UTl date) - 2451545.0 

(2.106) 

-*u — 

36525 


The actual equivalent expression which is coded is: 

h ^ =27r(I7Tl Julian day fraction) + 67310'. 54841 

+ 8640184'. 812866 T u + 0'.093104 Tj - 6' .2 x 10~ G T 3 (2.107) 


This expression produces a time, UT 1 , which tracks the Greenwich hour angle of the real Sun to 
within 16 m . However, it really is sidereal time, modified to fit our intuitive desire to have the Sun 
directly overhead at noon on the Greenwich meridian. Historically, differences of UT 1 from a uniform 
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measure of time, such as atomic time, have been used in specifying the orientation of the Earth. Note 
that this definition has buried in it the precession constant since it refers to the mean equinox of date. 

By the very definition of “mean of date” and “true of date”, nutation causes a difference in the 
hour angles of the mean equinox of date and the true equinox of date. This difference, called the 
“equation of equinoxes", is denoted by cue and is obtained as follows: 


where the vector 



(2.108) 


(2.109) 


is the unit vector, in true equatorial coordinates of date, toward the mean equinox of date. In mean 
equatorial coordinates of date, this same unit vector is just (1,0, 0) r . The matrix A r t -1 U the inverse 
(or equally, the transpose) of the transformation matrix N, which will be defined below in Eq. (2.120), 
to effect the transformation from true equatorial coordinates of date to mean equatorial coordinates 
of date. 

It is convenient to apply “ UXY ” as a group. To parts in 10 12 , XY = YX. However, with the 
same accuracy UXY / XYU . Neglecting terms of 0(© 2 ) (which produce station location errors of 
approximately 6 x 10 -4 cm): 


( cos H 
UXY = sin# 

\ sin ©i 

As for station coordinates, a time-linear model is also available for UTPM. If PM and UTl are 
symbolized by ©1-3, and the reference time is t 0 , then the model is 

©i = ©° + ©, (t - t 0 ) (2.111) 

where ©° are the values of UTPM at the reference epoch. 

2. 6. 1,1 Tidal UTPM Variations 

Tidal shifts of mass in the solid Earth, oceans, and atmosphere produce angular momenta 
which must be redistributed to satisfy conservation of total angular momentum. The consequences 
are variations in the orientation and rotation rate of the Earth: modification of polar motion and 
UTl. Such small effects emerged above the detection threshold in space geodesy in the early 1990 s. 
Modeling them is important if centimeter-level accuracy is to be obtained in interpretation of VLBI 
measurements. 

Just as various tidal forces affect the station locations (Secs. 2.4.1 - 2.4.4), they also affect polar 
motion and UTl (© 1 - 3 ). Equations similar to (2.45) may be written for each of the three components 
of Earth orientation: 

©1 = ©»o + A© lJO , + A ©, ocn -t- AQ tatm (2.112) 

where ©, (i— 1 , 3 ) symbolizes each of the three components of UTPM, © i0 is its value in the absence of 
tidal effects, and the three A terms are the respective contributions of solid Earth, ocean, and atmo- 
spheric tides. The next two sections describe the current models of solid and ocean tide contributions 
that are implemented in MODEST. At present, not enough is known about atmospheric tidal effects. 


— sin H — sin ©1 


( 2 . 110 ) 
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2.6.1. 1.1 Solid Earth Tide UTPM Variations 

The pioneering work in tidal effects on Earth orientation was that of Yoder et al. (1981). 
It was limited to UTl, but included some ocean effects. Their calculated AUTl can be represented 


as 


AUTl = 


N 


»=1 


Ai sin 


13 

■ 3 - 1 


(2.113) 


where N (=41) is chosen to include all terms with periods from 5 to 35 days. There are no other 
contributions until a period of 90 days is reached. However, these long-period terms are included by 
the measurements of the current Earth-orientation measurement services. The values for k+j and A*, 
along with the period involved, are given in Table VI. The for i = 1 , 5 are just the five fundamental 
arguments defined in Eqs. (2.123-2.127) as 2, F> D, and ft, respectively. In Table VI, the sign of 
the 14.73 day term has been changed [Yoder (1982)] to correct a sign error in Yoder et al. (1981). 
The BIH Annual Report for 1982 [BIH (1983)] is the first reference to give the correct table. 


2.6. 1.1.2 Ocean Tide UTPM Variations 


Redistribution of the angular momentum produced by ocean tides affects the Earth’s 
rotation pole position and velocity. This effect was first quantified by Brosche ct al. (1989, 1991). 
The dominant effects on polar motion and UTl are at diurnal, semidiurnal, fortnightly, monthly, 
and semiannual frequencies. Assuming that the frequencies slower than fortnightly are adequately 
accounted for either in modeling combined solid and ocean tidal effects (not strictly true with the 
Yoder model), or are already present in the a priori UTPM series, only the diurnal and semidiurnal 
frequencies need to be modeled. Further limiting the model to tidal components with apparent 
amplitudes larger than l*xs gives eight components. 

For unified notation, define ©i _ 3 = x, y polar motion and UT 1 , respectively. Then the ocean 
tidal effects A© can be written as 


N 

= V' 


A©< = 

*=i l 


r S 


An cos 


^3 kjOt} + + *) 

L J = 1 


Bn sin 


13 k 'i a ) 

L /=l 


ni(h 7 + r) 


(2.114) 


An and B t i are the cosine and sine amplitudes that may be calculated from theoretical tidal models 
(as in the work of Brosche) or determined from data (Herring and Dong, 1991; Herring, 1992; Sovers, 
Jacobs and Gross, 1993). Theoretical calculations of polar motion ocean effects have only very re- 
cently appeared: Gross (1993) used Seiler’s (1989-1991) ocean model to estimate 1 — 2 and Bi t i _2 

(i = 2 to 7 ). Table VII lists the eight terms currently included in the model. The numerical coefficients 
are taken from the results of Sovers ct al. (1993); this is known in MODEST code as the JPL92 fast 
UTPM model. 

The ocean tidal UTPM effects are also modulated by the 18. 6-year lunar node variation (N f ). 
As in the case of ocean loading station displacements (Sec. 2.4.3), the contributions A©/ of the 
companion tides to A 0 ^ can be written as 


N 


A0 ' , = EE rfc * 


i= 1 k 


An cos 


J] k t: a : + n t (h 7 + 1 r) + n*,*wjv*)t + n ki N‘ 


V* 1 


Bn sin 


T: k t yay + n,(/i^ + tt) + 


3 = 1 


(2.115) 


where the strengths of the companion tides r kt are found in Table V. These corrections are optionally 
available in MODEST. 

Since a rotational frequency of 1 ( 0 ) cycles per sidereal day (cpsd) in the celestial frame (5) 
is identical to a frequency of 0 (-1) cpsd in the Earth-fixed frame (B), nutations with space-fixed 
frequencies uj 5 coincide with polar motions with body-fixed frequencies — —1 + W5- The polar 
motion terms with coefficients A t 1,2 and B,i t 2 corresponding to the diurnal tidal components listed in 
Table VIII (1 = 5 to 8 ) are thus equivalent to components of the nutation model, and due care must 
be taken when both classes of parameters are estimated. 
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Table VI 

Periodic Tidally Induced Variations in UTl 
with Periods Less than 35 Days 


Index 

i 

Period 

(days) 

Argument coefficient 
fctl k\2 k»3 

A, 

(0V0001) 

1 

5.64 

1 

0 

2 

2 

2 

-0.02 

2 

6.85 

2 

0 

2 

0 

1 

-0.04 

3 

6.86 

2 

0 

2 

0 

2 

-0.10 

4 

7.09 

0 

0 

2 

2 

1 

-0.05 

5 

7.10 

0 

0 

2 

2 

2 

-0.12 

6 

9.11 

1 

0 

2 

0 

0 

-0.04 

7 

9.12 

1 

0 

2 

0 

1 

-0.41 

8 

9.13 

1 

0 

2 

0 

2 

-0.99 

9 

9.18 

3 

0 

0 

0 

0 

-0.02 

10 

9.54 

-1 

0 

2 

2 

1 

-0.08 

11 

9.56 

-1 

0 

2 

2 

2 

-0.20 

12 

9.61 

1 

0 

0 

2 

0 

-0.08 

13 

12.81 

2 

0 

2 

-2 

2 

0.02 

14 

13.17 

0 

1 

2 

0 

2 

0.03 

15 

13.61 

0 

0 

2 

0 

0 

-0.30 

16 

13.63 

0 

0 

2 

0 

1 

-3.21 

17 

13.66 

0 

0 

2 

0 

2 

-7.76 

I 18 

13.75 

2 

0 

0 

0 

-1 

0.02 

19 

13.78 

2 

0 

0 

0 

0 

-0.34 

20 

13.81 

2 

0 

0 

0 

1 

0.02 

21 

14.19 

0 

-1 

2 

0 

2 

-0.02 

22 

14.73 

0 

0 

0 

2 

-1 

0.05 

23 

14.77 

0 

0 

0 

2 

0 

-0.73 

24 

14.80 

0 

0 

0 

2 

1 

-0.05 

25 

15.39 

0 

-1 

0 

2 

0 

-0.05 

26 

23.86 

1 

0 

2 

-2 

1 

0.05 

27 

23.94 

1 

0 

2 

-2 

2 

0.10 

28 

25.62 

1 

1 

0 

0 

0 

0.04 

29 

26.88 

-1 

0 

2 

0 

0 

0.05 

30 

26.98 

-1 

0 

2 

0 

1 

0.18 

31 

27.09 

-1 

0 

2 

0 

2 

0.44 

32 

27.44 

1 

0 

0 

0 

-1 

0.53 

33 

27.56 

1 

0 

0 

0 

0 

-8.26 

34 

27.67 

1 

0 

0 

0 

1 

0.54 

35 

29.53 

0 

0 

0 

1 

0 

0.05 

36 

29.80 

1 

-1 

0 

0 

0 

-0.06 

37 

31.66 

-1 

0 

0 

2 

-1 

0.12 

38 

31.81 

-1 

0 

0 

2 

0 

-1.82 

39 

31.96 

-1 

0 

0 

2 

1 

0.13 

40 

32.61 

1 

0 

•2 

2 

-1 

0.02 

41 

34.85 

-1 

-1 

0 

2 

0 

-0.09 
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Table VII 

Ocean Tidally Induced Periodic Variations in Polar Motion (JPL92 Model) 


Index 

i 

Tide 

period 

(hours) 

kii 

Argument coefficient 
k\2 ^t3 

n : 

An 

Bn A{ 2 

(fias) 

Bi2 

1 

k 7 

11.967 

0 

0 

0 

0 

0 

-2 

-2 

65 

44 

-57 

2 

S 2 

12.000 

0 

0 

2 

-2 

2 

-2 

101 

166 

126 

-89 

3 

m 2 

12.421 

0 

0 

2 

0 

2 

-2 

26 

283 

247 

-2 

4 

n 2 

12.658 

1 

0 

2 

0 

2 

-2 

-15 

56 

19 

-11 

5 

K i 

23.934 

0 

0 

0 

0 

0 

-1 

-583 

-2780 

-2950 

376 

6 

P i 

24.066 

0 

0 

2 

-2 

2 

-1 

154 

46 

42 

-17 

7 

O ! 

25.819 

0 

0 

2 

0 

2 

-1 

242 

-152 

2 

-30 

8 

Qi 

26.868 

1 

0 

2 

0 

2 

-1 

72 

-32 

26 

7 


Table VIII 

Ocean Tidally Induced Periodic Variations in UTl (JPL92 Model) 


Index 

i 

Tide 

period 

(hours) 

kn 

Argument coefficient 
&t2 &i3 ^»4 


Aiz Biz 

(O.l/is) 

1 

k 2 

11.967 

0 

0 

0 

0 

0 

-2 

-9 

26 

2 

S 2 

12.000 

0 

0 

2 

-2 

2 

-2 

-4 

52 

3 

m 2 

12.421 

0 

0 

2 

0 

2 

-2 

-104 

149 

4 

n 2 

12.658 

1 

0 

2 

0 

2 

-2 

-23 

20 

5 

K x 

23.934 

0 

0 

0 

0 

0 

-1 

35 

151 

6 

Pi 

24.066 

0 

0 

2 

-2 

2 

-1 

-32 

-64 

7 

Oi 

25.819 

0 

0 

2 

0 

2 

-1 

-135 

-166 

8 

Q i 

26.868 

1 

0 

2 

0 

2 

-1 

-40 

-53 


2.6. 1.1. 3 Interpolation of UTPM Values 


Depending on the smoothing used to produce the a priori UTl — UTC series, the short- 
period (t < 35 days) fluctuations in UTl due to changes in the latitude and size of the mean tidal 
bulge may or may not be smoothed out. Since we want as accurate an a priori as possible, it may be 
necessary to add this effect to the UTl o priori obtained from the series UTl 9rnoo thtd • R this option 
is selected, then the desired a priori UTl is given by 

UTla priori = UTl.moothtd + AC/Tl (2.116) 

UTl 9 moothed represents an appropriately smoothed a priori measurement of the orientation of the 
Earth (i.e., typically IERS Bulletin A smoothed or, even better, UTlR), for which the short period 
(t < 35 days) tidal effects have either been averaged to zero, or, as in the case of UTlR } removed 
before smoothing. 

It might be appropriate at this point to describe the interpolation method used in MODEST 
to obtain a prion polar motion and UTl values. These are normally available as tables at 5-day or 
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1-day intervals, from either IERS ( e.g . IERS, 1993) or the International Radio Interferometric Sur- 
veying (IRIS) project of the International Association for Geodesy (IAG) ( t.g . IAG, 1993). Linear 
interpolation is performed for all three quantities. If the short-period tidal terms AITTl are present 
in the tabular values, they are subtracted before interpolation, and added back to the final value. 
With the present accuracy of determinations of pole position and UT 1 (<1 mas and 0.05 ms respec- 
tively), linear interpolation over a 5-day interval may be inadequate, possibly giving rise to 0.1 ms 
errors in UTl. Quadratic spline interpolation is a possible alternative. Even with the present code, 
however, the highest possible accuracy may be achieved by performing the interpolation externally 
to MODEST, and supplying it with tables of values more closely spaced in time for the final internal 
linear interpolation. The Kalman-filtered UTPM values of Gross (1992), with values given at 1-day 
intervals, are ideally suited for this purpose. 
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2.6.2 NUTATION 


With the completion of the UT 1 and polar motion transformations, we are left with a station 
location vector, This is the station location relative to true equatorial celestial coordinates of 

date. The last set of transformations are nutation, N , precession, P, and the perturbation rotation, 
H, applied in that order. These transformations give the station location, r c , in celestial equatorial 
coordinates: 

Tc-nWVr**. (2.117) 

The transformation matrix TV is a composite of three separate rotations (Melbourne tt a/,, 1968): 

1. A(c): true equatorial coordinates of date to ecliptic coordinates of date. 

(10 0 \ 

A(e) = 0 cos c sine I (2.118) 

\0 — sin c cos e J 

2. C T (Sip): nutation in longitude from ecliptic coordinates of date to mean ecliptic coordinates of 
date. 

( cos Sip sin Sip 0\ 

— sinSip cos Sip 0 j (2.119) 

0 0 ! J 

where Sip is the nutation in ecliptic longitude. 

3. A T (?): ecliptic coordinates of date to mean equatorial coordinates. 


In ecliptic coordinates of date, the mean equinox is at an angle Sip = tan” 1 (t/^/x^)- The angle 
Se = c — ? is the nutation in obliquity, and ? is the mean obliquity (the dihedral angle between the 
plane of the ecliptic and the mean plane of the equator). “Mean* as used in this section implies that 
the short-period ( T < 18.6 years) effects of nutation have been removed. Actually, the separation 
between nutation and precession is rather arbitrary, but historical. The composite rotation is: 

TV = A T (e)C T {6iP)A{c) (2.120) 


( cos Sip 
— cos ?sin Sip 
— sin ? sin Sip 


cos e sin Sip 

cos? cos e cos Sip 4- sin? sine 
sin ? cos e cos <5 ip — cos? sine 


sin esin Sip 
cos? sine cos Sip — sin? cose 
sin ? sin e cos Sip + cose cose 


The 1980 IAU nutation model (Seidelmann, 1982; Kaplan, 1981) is used to obtain the values for 
Sip and e — ?. The mean obliquity is obtained from Lieske ct aL (1977) or from Kaplan (1981): 


e = 23° 26' 2l/ ; 448 - 46. "8150 T - 5. "9 x 10 _4 T 2 + l"813 x 10" 3 T 3 (2.121) 


_ (Julian TDB date) — 2451545.0 
” 36525 


(2.122) 


This nutation in longitude (Sip) and in obliquity ( Se = e — ? ) can be represented by a series expansion 
of the sines and cosines of linear combinations of five fundamental arguments. These are (Kaplan, 
1981; Cannon, 1981): 

1. the mean anomaly of the Moon: 


Qi = / = 485866". 733 + (I325 r + 715922".633) T 

+ 31". 310 T 2 + 0".064 T 3 (2.123) 

2. the mean anomaly of the Sun: 

an = V = 1287099". 804 + (99 r + 129258l".224) T 

- 0".577 T 2 - 0."012 T 3 (2.124) 
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3. the mean argument of latitude of the Moon: 

q 3 = F = 335778". 877 + (1342 r + 295263".137) T 

- 13". 257 T 7 + 0".011 T 3 (2.125) 

4 . the mean elongation of the Moon from the Sun: 

q 4 = D = 1072261". 307 + (1236 r + 1105601".328) T 

- 6". 891 T 2 + 0".019 T 3 (2.126) 

5 . the mean longitude of the ascending lunar node: 

a 5 = n = 450160". 280 - (5 r + 482890".539) T 

+ 7". 455 T 2 + 0".008 T 3 (2.127) 


where l r = 360° = 1296000". 

With these fundamental arguments, the nutation quantities can then be represented by 


N 

r * 

* 

w-E 

(A 0j + AijT) sin y^fc^ai(T) 


y«x 

. isl 

- 


N 

* 

r 5 1 

m 


(Boy + B\ yT) cos 

E^ a -( r ) 


y=i 

. 

s=i 1 

. 


where the various values of a^, kji , Aj t and By are tabulated in Table A. I. 
2.6. 2.1 Corrections to the 1980 IAU Model 


(2.128) 


(2.129) 


Additional terms can be optionally added to the nutations Sip and 6c in Eqs. (2.128) 
and (2.129). These include the out-of-phase nutations, the free-core nutations (Yoder, 1983) with 
period c oj (nominally 430 days), and the “nutation tweaks” Aip and Ac, which are arbitrary constant 
increments of the nutation angles 6rp and 6c. Unlike the usual nutation expressions, the tweaks have 
no time dependence. The out-of-phase nutations, which are not included in the IAU 1980 nutation 
series, are identical to Eqs. (2.128) and (2.129), with the replacements sin cos: 


and 


N 

y=i 


• 

5 -j‘ 

{A 2} + A Z jT ) cos 

k=i ■». 


N 


^=E 

j=i 


(-Bay + A/T) sin 


E m 


*■ 1=1 


(2.130) 

(2.131) 


Expressions similar to these are adopted for the free-core nutations: 

S^ } - (/loo + .4io7') sin(w;T) + {A 20 + A 20 T) cos( W/ r) (2.132) 

and 

6e^ = (B 0 o + BioT) cos (u?jT) + (B 20 + B 30 T) sin(a>/T) (2.133) 

If the free-core nutation is to be retrograde, as expected on theoretical grounds, uj should be negative. 
The nutation model thus contains a total of 856 parameters: A t0 (i=0,3; ; = 1,106) and B i; (i =0,3 ; 
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j= 1,106) plus the free-nutation amplitudes A^o (*=0,3), Bio (t=0,3). The only nonzero a priori 
amplitudes are the A 0 y, A ljt Boy, Biy (7 = 1,106) given in Table A.I. 

The nutation tweaks are just constant additive factors to the angles Sip and 6c: 

Sip —* Sip + Aip (2.134) 


fie — * + Ac (2.135) 


Several alternatives are available as MODEST options to correct deficiencies in the IAU nutation 
model. The first possibility is to use empirically determined values of A ip, A c as part of the polar 
motion and UTl input which was described in Section 2.6.1. 1.3. If this option is selected, the user is 
relying on nutation angles that are determined from other VLBI experiments near the date of interest, 
and performing linear interpolation. 

Other options are available to select one of the recently published replacements of the 1980 IAU 
series. Zhu et al. (1989, 1990) have refined the 1980 IAU theory of nutation both by reexamining the 
underlying Earth model and by incorporating recent experimental results. The Zhu et al. results are 
tabulated here in three parts: a) the original 106 terms of the 1980 IAU series with revised amplitudes 
in Table A. II, b) four sets of out of phase terms in Table A.III, and c) an additional 156 terms in 
Table A.IV. Herring (1991) has extended the work of Zhu ct al. and used geophysical parameters from 
Mathews et al. (1991) to generate the ZMOA 1990-2 (Zhu, Mathews, Oceans, Anelasticity) nutation 
series. Tables A.V to A. VII list the coefficients of this series. Kinoshita and Souchay (1990) have 
reexamined the rigid-Earth nutation theory, and calculated all terms larger than 0.005 mas, including 
planetary terms not included in any previous theories. The 263 lunisolar terms listed in Table A. VIII 
have been corrected for Earth non-rigidity (Souchay, 1993); Table A.IX lists the 117 planetary terms. 

Note that the paper of Kinoshita and Souchay gives expressions for the lunisolar tidal arguments 
that are at variance with the IAU formulas presented above in Eqs. (2.123-2.127). Their expressions 
for the five usual arguments l , l 1 , F , D, and 0, as well as five additional planetary arguments lv , Ibi 
l M , h, and Is are given below (all in units of radians): 

1. the mean anomaly of the Moon: 

ai = / = 2.35555590 + 8328.691427 T (2.136) 

2. the mean anomaly of the Sun: 

q 2 = /' = 6.24006013 + 628.301955 T (2.137) 

3. the mean argument of latitude of the Moon: 

a 3 = F= 1.62790523 + 8433.466158 T (2.138) 

4. the mean elongation of the Moon from the Sun: 

a 4 = D = 5.19846674 + 7771.377147 T (2.139) 

5. the mean longitude of the ascending lunar node: 

= fl = 2.18243920 - 33.757045 T (2.140) 


The Kinoshita-Souchay planetary contributions to Sip and Se are 

n r rio 1 io 

Sip = 7: S^j sin X>y,&m + C^j cos Ekyiftm 

; = 1 L L »= I L »=l 
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and 


N 

w 

r 10 -I 

r 10 1 

■ 


Sgj sin 


+ cos ly^ftift(T) 


/= 1 

. 

1 

k=i 

. 


where the astronomical arguments are symbolized by ft; the last four ft are identical with the as 
defined above (ft = D — a 4 , ft = F = 03, ft = / = ai, fto = O = as), while the first six are 

1. the mean anomaly of Venus: 


ft = lv = 3.176146697 H- 1021.3285546 T (2.143) 

2. the mean anomaly of Earth: 

ft = i E = 1.753470314 + 628.30758492 T (2.144) 

3. the mean anomaly of Mars: 

ft = l M = 6.203480913 + 334.06124315 T (2.145) 

4. the mean anomaly of Jupiter: 

ft lj = 0.599546497 + 52.96909651 T (2.146) 

5. the mean anomaly of Saturn: 

ft = l s = 0.874016757 + 21.32990954 T (2.147) 

6. general precession: 

ft = Pa = 0.02438175 T + 5.38691 x 10‘ 6 T 2 (2.148) 

For simulation purposes, the older Woolard nutation model is also available in MODEST. With 
the exception of the number, amplitudes, and arguments of the terms, the older series is exactly 
analogous to the 1980 IAU theory, i.e., of the form of Eqs. (2.128) and (2.129). For completeness of 
documentation, the coefficients are listed in Table A.X. 

No partial derivatives with respect to the Woolard or Zhu et al. amplitudes are currently calcu- 
lated. It is emphasized that, for the present, the default nutation model in MODEST remains as just 
the 1980 IAU nutation model given in Table A. I. 
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2.6.3 PRECESSION 


The next transformation in going from the terrestrial frame to the celestial frame is the rotation 
P. This is the precession transformation from mean equatorial coordinates of date to the equatorial 
coordinates of the reference epoch (e.g. t J2000). It is a composite of three rotations discussed in detail 
by Melbourne ct at. (1968) and Lieske et al. (1977): 


R(-Z) = 


cos Z sin Z 0 
— sin Z cos Z 0 
0 0 1 


(2.149) 


( cos © 0 sin © 

0 1 0 | (2.150) 

— sin © 0 cos © 

( cos £ sin £ 0 \ 

— sin £ cos £ 0 (2.151) 

0 0 l) 

P = R(-s)Q(e)R{-Z) (2.152) 

( cos f cos © cos Z - sin £ sin Z cos £ cos ©sin Z 4 sin £ cos Z cos £ sin © 

— sin £ cos © cos Z — cos £ sin Z — sin £ cos © sin Z 4 cos f cos Z — sin $ sin © 

— sin © cos Z — sin © sin Z cos © 

The auxiliary angles ©, Z depend on precession constants, obliquity, and time as 

f = 0".5mT 4 0". 30188 T 2 4 0".017998 T 3 (2.153) 

Z = 0".5mT 4- l". 09468 T 2 4- 0".018203 T 3 (2.154) 

© = nT - 0". 42665 T 2 - 0".041833 T 3 (2.155) 

where the speeds of precession in right ascension and declination are, respectively, 

m *= Pls cos £o ~ Ppl (2.156) 

n = pls sin £o (2.157) 


and pls = the luni-solar precession constant, ppl = planetary precession constant, £o — the obliquity 
at J2000, and T [Eq. (2.122)] is the time in centuries past J2000. Nominal values at J2000 are pls 
= 5038".7784/cy, ppl — 10".5526/cy; these yield the expressions given by Lieske et al. (1977) and 
Kaplan (1981): 


$ = 2306".218lT 4 0". 30188 T 2 4 0".017998 T 3 (2.158) 

© = 2004". 3109 T - 0". 42665 T 2 - 0".041833 T 3 (2.159) 

Z = 2306 ,/ .2181 T 4 l ,/ .09468 T 2 4 0".018203 T 3 (2.160) 

Partial derivatives of the VLBI observables with respect to luni-solar and planetary precession are 

derived from the expressions (2.152-2,157) and given in Section 2.9. The precession matrix completes 
the standard model for the orientation of the Earth. Numerical checks of direct estimates of preces- 
sion corrections against similar estimates based on the perturbation rotation (next section) ensure 
consistency. 

2.6.4 PERTURBATION ROTATION 

This standard model for the rotation of the Earth as a whole may need a small incremental 
rotation about any one of the resulting axes. Define this perturbation rotation matrix as 

n = A T A y A 2 (2.161) 
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where 


A, - 



(2.162) 


with 6G, being a small angle rotation about the x axis, in the sense of carrying y into z; 

1 0 -se v 


A„ = 


se v o 


0 

1 


(2.163) 


with 6G V being a small angle rotation about the y axis, in the sense of carrying z into x; and 


f 1 6G, 0 

A, = -SG, 1 0 

V 0 0 1 


(2.164) 


with 6G, being a small angle rotation about the z axis, in the sense of carrying x into y. For angles 
on the order of 1 arc second we can neglect terms on the order 6G 7 R E as they give effects on the 
order of 0.015 cm. Thus, in that approximation 

f 1 se, —se v \ 

n = -se, i se, (2.105) 

V se y -se, 1 J 

In general, 

SGi = £ 6 <(t) = 6 e i0 + SGiT + fi[T) (2.166) 

which is the sum of an offset, a time-linear rate, and some higher order or oscillatory terms. Currently, 
only the offset and linear rate are implemented. In particular, a non-zero value of 6G V is equivalent to 

a change in the precession constant, and 6G, is equivalent to the time rate of change of the obliquity 
e. Setting 

se, = se v = se, = o (2.167) 

gives the effect of applying only the standard rotation matrices. 

Starting with the Earth-fixed vector, To, we have in Sections 2.3 through 2.6 above shown how 
we obtain the same vector, r c , expressed in the celestial frame: 


t c = {1PNUXY(tq + A) 


(2.168) 


2.7 EARTH ORBITAL MOTION 

We now wish to transform these station locations from a geocentric reference frame moving 
with the Earth to a celestial reference frame which is at rest relative to the center of mass of the 
Solar System. In this Solar System barycentric (SSB) frame we will use these station locations to 
calculate the geometric delay (see Section 2.1). We will transform the resulting time interval back to 
the frame in which the time delay is actually measured by the interferometer - the frame moving 
with the Earth. 

Let E' be a geocentric frame moving with vector velocity = /3c relative to a frame, E, at rest 
relative to the Solar System center of mass. Further, let r(t) be the position of a point ( e.g ., station 
location) in space as a function of time, t, as measured in the E (SSB) frame. In the E' (geocentric) 
frame, there is a corresponding position r'(f') as a function of time, t'. We normally observe and 
model r'(f') as shown in Sections 2.3 through 2.6. However, in order to calculate the geometric delay 
in the SSB frame (E), we will need the transformations of r(t) and r'(t'), as well as of t and t', as 
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we shift frames of reference. 
(1975): 


Measuring positions in units of light travel time, we have from Jackson 

r'(t') = r(t) + ft - l)r(t) • - ifit (2.169) 

P 

f = -r[t - r(t) ■ fi\ (2.170) 


and for the inverse transformation: 




(2.171) 


t = -r[t'+r'(0-/9] 


(2.172) 


where 


l = {l-fi 7 ) 


- 1/2 


(2.173) 


Let ti represent the time measured in the SSB frame (E), at which a wave front crosses antenna 
1 at position r^tj). Let r 2 (ti) be the position of antenna 2 at this same time, as measured in the 
SSB frame. Also, let t‘ 2 be the time measured in this frame at which that same wave front intersects 
station 2. This occurs at the position r 2 (^ 2 ). Following Section 2.1, we can calculate the geometric 
delay t 2 — ti- Transforming this time interval back to the geocentric (E') frame, we obtain 


f 2 ' - ti = Tf(t; - ti) - n\T,(ti) - Ti (ti)] • fi 


(2.174) 


Assume further that the motion of station #2 is rectilinear over this time interval. This assumption 
is not strictly true but, as discussed below, the resulting error is much less than 1 mm in calculated 
delay. Thus, 

r 2(*5) = r 2(*i) + 02^2 * 1 ) (2.175) 


which gives: 

r 2 (t 2 ) - ri(ti) = r 2 (ti) - ri(*i) +£ 2 (t£ - ti) 

and 


(2.176) 


t m 2 - t\ = l{t* 2 - h) - l\ri(t x ) - ri(ti)] 0-70 2 ’ fi \ l 2 ~ *i] 

= 7(1 - 02 • fi)i$ - *1) - - n(*i)] ■ 0 ( 2 - 177 ) 


This is the expression for the geometric delay that would be observed in the geocentric (E') frame in 
terms of the geometric delay and station positions measured in the SSB system (E). 

Since our calculation starts with station locations given in the geocentric frame, it is convenient 
to obtain an expression for [r 2 (ti) — r x ( ^ 1 ) ] terms of quantities expressed in the geocentric frame. 

To obtain such an expression consider two events [r x (^i ) * r 2(*i)l that are geometrically separate, but 
simultaneous, in the geocentric frame, and occurring at time t[ . These two events appear in the SSB 
frame as: « « 

ri(t») = r' (*' ) + (7 - l)ri(ti) ■ y + lfit\ (2-178) 

and as: «« 

r 2 (h) = r' 7 (t\) + ft - l)r 2 (t'j) ■ % + ifit\ (2.179) 

P 

where 

* a -ti=7[ri(t;)-ri(t , 1 )]^ (2.180) 

With these three equations and the expression 

r2 (t 2) = 1*0 (t 1 ) + 02^2 — *i] ( 2 . 181 ) 
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we may obtain the vector r 2 (ti): 


r a (M = ri(ti) + (y ~ IWW) • S’ + 101 ~ iMrWi) " P (2-182) 

P 

This is the position of station #2 at the time t x as observed in E. PVom this we obtain: 

r*(ti) - rx( tl ) = r' (t'J - r '(«') + (7 - l)|r' a (*i) - r'(t')] ^ 

-^ 2 K(£'x)-r' x (ti)] P (2-183) 

As shown in Section 2.1, the vectors (r 2 (*i) — ri (ti)] and 0 2 are all that is needed to obtain t\ - t x for 
the case of plane waves. For curved wave fronts we will need to know the individual station locations 
in the barycentric frame as well. These we obtain from (2.178) and (2.182) with t[ set equal to zero. 
Setting ti s= 0 is justified since the origin of time is arbitrary when we are trying to obtain time 
differences. 

In the actual coding of these transformations, the relationship for the transformation of velocities 
is also needed. Taking differentials of (2.171) and (2.172) we have; 


dr = dr 1 4 (7 — l)dr' * 4- ~iPdt* 

P 

dt = *i(dt’ 4 dr' * P) 

Dividing to obtain dr/dt we obtain for station #2 in the SSB (E) frame: 


A + b-iJA-y + 7* 

^ 2= 7(1 +P' 7 P) 

For station #2 relative to the geocentric origin, we have from (2.98) and (2.99): 

r dU 


where 


P \ « n pn—xyt^wb 


w B = 7.2921151467 x 10" 5 rad/s 


(2.184) 

(2.185) 

(2.186) 

(2.187) 

(2.188) 


is the inertial rotation rate of the Earth as specified in Kaplan (1981), p.12. This is not a critical 
number since it is used only for station velocities, or to extrapolate Earth rotation forward for very 
small fractions of a day (x.e., typically less than 1000 seconds). Actually, this expression is a better 

dH 

approximation than it might seem from the form since the errors in the approximation, = W£, 

are very nearly offset by the effect of ignoring the time dependence of PN . 

The assumption of rectilinear motion can be shown to result in negligible errors. Using the plane 
wave front approximation (2.2), we can estimate the error St in the calculated delay due to an error 
A P 2 in the above value of P 2 : 


St = k [r 2 (ti) - ri(ti) 
Further, from (2.186) above, 


since 


1 — k (f ) 2 4 A^ 2 ) 1 “ h fi 2 

AP 2 « &p 2 

^14 10 -8 


£3 T A^ 2 


(2.189) 

(2.190) 

(2.191) 
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For the vector in a frame rotating with angular velocity w, the error A#* that accumulates in the 
time interval r due to neglecting the rotation of that frame is 

A^2 « ^2 u;r (2.192) 


Thus for typical Earth-fixed baselines, where r < 0.02 s, neglect of the curvilinear motion of station 
#2 due to the rotation of the Earth causes an error of < 4 x 10” 14 s, or 0.0012 cm, in the calculation 
of r. Similarly, neglect of the orbital character of the Earth’s motion causes an error on the order of 
0.00024 cm maximum. 

The position, Re, and velocity, 0 E , of the Earth’s center about the center of mass of the Solar 
System are: 

Re = (2.193) 


& E 




(2.194) 


where the index i indicates the Sun, Moon, and all nine Solar System planets, m, is the mass of the 
body indexed by t, while R, and are that body’s center-of-mass position and velocity relative to 
the center of the Earth in the barycentric frame. In a strict sense, the summation should be over 
all objects in the Solar System. Except for the Earth-Moon system, each planet mass represents not 
only that planet’s mass, but also that of all its satellites. The Ri and ^ are obtained from the JPL 
planetary ephemeris (DE200 as of May, 1982) for the J2000 frame. 

Working in a frame at rest with respect to the center of mass of the Solar System causes rel- 
ativistic effects due to the motion of the Solar System in a “fixed frame” to be included in the 
mean position of the sources and in their proper motion. The effects of galactic rotation can 
be easily estimated. In the vicinity of the Sun, the period for galactic rotation is approximately 
2.2 xlO 8 years. Thus our angular velocity about the galactic center is « 2?r/2.2 x 10® = 3 x 10 8 ra- 
dians/year. For sources within the Galaxy, at distances approximately equal to our distance from the 
galactic center, therefore, the apparent positions could change by « 30 nrad/yr. An intercontinental 
baseline (10,000 km) could thus be in error by as much as 30 cm/yr (1 nrad « 1 cm) if measurements 
were based on sources within the Galaxy. Since our distance from the galactic center is « 2.7 x 10 4 
light years, and most extragalactic radio sources are believed to be « 10 9 light years distant, the 
potential baseline error is scaled by the ratio of these two distances, «3x 10” 5 , and becomes « 0.001 
cm/year. Even with the present 15-year history of VLBI data, the purely geometric systematic error 
due to galactic rotation is probably negligible, and only exceeds the millimeter level for sources closer 
than 100 million light years. 

A second contribution to time variability is aberration due to special relativistic effects for ob- 
servers on a moving platform. Both the galactic latitude and longitude of the Solar System vary 
sinusoidally with the galactic rotation period T = 2.2 x 10 8 years, amplitudes ranging over ±vfe 
(French, 1968), and the latitude variation is also proportional to the sine of the latitude. The ampli- 
tude range is 7.5 x 10~ 4 ; thus over half the rotational period at zero galactic latitude, the longitude 
varies by 

AB = (2u/c)/(0.5T) = 1.4 x 10" u rad/y (2.195) 


Over the present 15-year span of VLBI data, the systematic errors induced by aberration are thus 
0.2 nrad (40 /ias) in angular measurements, and 2 mm in distance measurements for a 10,000 km 
baseline. Both are in the range that is currently starting to be detectable, and serious consideration 
will be given to model enhancements in the near future. 
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2.8 ANTENNA GEOMETRY 


The above work indicates how the time delay model would be calculated for two points fixed 
with respect to the Earth’s crust. In practice, however, an antenna system does not behave as an 
Earth-fixed point. Not only are there instrumental delays in the system, but portions of the antenna 
move relative to the Earth. To the extent that instrumental delays are independent of the antenna 
orientation, they are indistinguishable to the interferometer from clock offsets and secular changes 
in these offsets. If necessary, these instrumental delays can be separated from clock properties by 
a careful calibration of each antenna system. That is a separate problem, treated as a calibration 
correction (e.g., Thomas, 1981), and will not be addressed here. 

However, the motions of the antennas relative to the Earth’s surface must be considered since 
they are part of the geometric model. A fairly general antenna pointing system is shown schematically 
in Figure 5. The unit vector, 8, to the apparent source position is shown. Usually, a symmetry axis 
AD will point parallel to s. The point A on the figure also represents the end view of an axis which 
allows rotation in the plane perpendicular to that axis. This axis is offset by some distance H from a 
second rotation axis BE. All points on this second rotation axis are fixed relative to the Earth. 
Consequently, any point along that axis is a candidate for the fiducial point which terminates this 
end of the baseline. The point we actually use is the point P. A plane containing axis A and perpen- 
dicular to BE intersects BE at the point P. This is somewhat an arbitrary choice, one of conceptual 
convenience. 



Figure 5. A generalized schematic representation of the geometry of a steerable antenna 


Consider the plane Q which is perpendicular to the antenna symmetry axis, AD, and contains 
the antenna rotation axis A. For plane wave fronts this is an isophase plane (it coincides with the 
wave front). For curved wave fronts this deviates from an isophase surface by sr H 2 /(2R), where R 
is the distance to the source, and H is taken as a typical antenna offset AP. For H ^ 10 meters, 
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D _ p_ _ 60 Rf « 3 6 x 10 8 m, and the curvature correction U 2 /(2R) « 14 x 10 J m is 
totally neltUeible R has to be 5 km, or 10~ 3 R E , before this deviation approaches 1 cm contribution. 
Consequently, for all anticipated applications of radio interferometry using high-gam ^d.o antennas 
the curvature of the wave front may be neglected in obtaining the effect on the time delay of the 

^^P^ov^dirthTkistrumental delay of the antenna system is independent of the antenna orientation, 
the recorded signal is at a constant phase delay, independent of antenna orientation, at any poin on 
the Q plane. Since this delay is indistinguishable from a clock offset, it will be totally absorbed by 

that portion of our model. 

2.8.1 AXIS OFFSET 

The advantage of choosing the Q plane rather than some other plane parallel to it is that axis 
A is contained in this plane, and axis A is fixed relative to the BE axis by the antenna structure^ If I 
is the length of a line from P perpendicular to the Q plane, the wave front will reach the Earth-fixed 
point P at a time At = l/e after the wave front passes through axis A. If r 0 is ^e model delay for a 
wave front to pass from P on antenna #1 to a similarly defined point on antenna #2, then the model 
for the observed delay should be amended as: 

T = r 0 - (At 2 - Ati) = To + (*i - k)/e (2.196) 


where the subscripts refer to antennas #1 and #2. ...... ,, Q n c . 

For the inclusion of this effect in the model, we follow a treatment given by Wade (1970). Define 
a unit vector I along BE, in the sense of positive away from the Earth. Further, define a vector, L, 
from P to A. Without much loss of generality in this antenna system, we assume that s,L, and I are 

coplanar. Then: 


L = ±# £ 


Ixjsxl] 


Ix[s X I] I 


(2.197) 


where the plus or minus sign is chosen to give L the direction from P to A. The plus sign is used if, 
when 8 and L are parallel or antiparallel, the antenna comes closer to the source as H increases. Since 

(2.198) 


Ixjsxl] = 8- I [i s] 

l = b ■ L = ±H yl — je • I] 


(2.199) 


where the sign choice above is carried through. ^ , . , _ . 

Curvature is always a negligible effect in the determination of e L. Likewise, gravitational effects 
are sufficiently constant over a dimension jL| so as to enable one to obtain a single Cartesian frame 
over these dimensions, to a very good approximation. Consequently, it is somewhat easier to calculate 
a proper time A t-l/c in the antenna frame and to include it in the model by adding it to r 0 , taking 
into account, in principle at least, the time dilation in going from the antenna frame to the frame in 
which tq is obtained. 


2.8.2 REFRACTION 

Thus, if So is the unit vector to the source from the antenna in a frame at rest with respect to 
the Solar System center of mass, perform a Lorentz transformation to obtain s, the apparent source 
unit vector in the Earth-fixed celestial frame. Actually, the antenna does not “look” at the apparent 
source position s, but rather at the position of the source after the ray path has been refracted 
by an angle e in the Earth’s atmosphere. This effect is already included in the tropospheric delay 
correction (Section 4); however since the antenna model uses the antenna elevation angle E 0 , the 
correction must be made here as well. For the worst case (elevation angle of 6°) at average DSN 
station altitudes, the deflection can be as large as 2 x 10 -3 radians. Thus, 61 « He as 2 cm for 
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^ = ^0 meters. A model option permits modification of so to take atmospheric refraction into 
account. The large-elevation-angle approximation is the inverse tangent law: 

A E = 3.13 x 10“ 4 /tanjE o (2.200) 

where E is the elevation angle, and A E the change in apparent elevation E 0 induced by refraction. This 
model was implemented only for software comparison purposes, since it gives incorrect results at low 
elevation angles. In the notation of Section 4.2, a single homogeneous spherical layer approximation 
yields the bending correction in terms of the zenith troposphere delays pz> refractivity moment A/nni 


scale height A, and Earth radius R : 

L.E = cos _1 [cos(£ 0 + c*o)/(l + xo)] - e» 0 (2.201) 

where 

Xo = {pz itv + Pz w „/M 0 oi)/A (2.202) 

a 0 = cos -1 [(l + cr')/(l + <r)] (2.203) 

<r = A/R (2.204) 

o' = [ (l + o(a + 2)/ sin 2 Eo) — lj sin 2 Eq (2.205) 


This formula agrees with ray-tracing results to within 1% at 6° and esl5% at 1° elevation, while the 
corresponding comparisons for Eq. (2.200) give «25% at 6° and a factor of 3 at 1°. 

Since we are given I in terrestrial coordinates, we first perform the coordinate transformation 
given by Q above: 

I = QIterr««tr*o/ (2.206) 

With this done, obtain At = l/e , as shown in Figure 6 for each of the major antenna types. Note 
that for “nearby" sources we also must include parallax (e.g., geographically separate antennas are 
not pointing in the same direction). If Ro is the position of the source as seen from the center of the 
Earth, and r is the position of a station in the same frame, then the position of the source relative to 
that station is 

R = Ro - r (2.207) 

and in (2.199) we make the substitution 



‘ [Ro-r] l 1 2 

. |Ro-r| _ 


(2.208) 


2.8.S UNIQUE ANTENNAS 


One of the VLBI antennas employed by the IRIS project of the National Geodetic Survey 
does not fall into any standard category. It is unique because it is an equatorial mount designed for 
the latitude of Washington, D.C., but was deployed at Richmond, Florida until it was destroyed in the 
hurricane of August 1992. The considerable latitude difference, and the axis offset of several meters, 
make it imperative that the antenna geometry be properly modeled. In the local VEN coordinate 
frame, the vector I is 

( sin 4> w A 

I - cos <t> w sine (2.209) 

V COS 4>w cos € J 


Upon transformation to the Earth-fixed frame via the matrix VW (Eq. (2.71)], it becomes 

( cos A (sin 4>w cos 4> — cos <f>w sin cos c) + sin A cos <f>w sin e \ 
sin A(sin <pw cos 4> - cos 4>w sin 4> cos e) - cos A cos <f>w sin c 
sin <f>w sin 4> -f- cos 4>w cos <f> cos e J 


( 2 . 210 ) 
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Here (A,^) are the Richmond longitude and latitude, <j>w is the latitude of Washington (39.06°), and 
€ = 0.12° W of N is the azimuth misalignment. 

Two other one-of-a-kind antennas, Arecibo and Nancay, are seldom used in astrometric and 
geodetic VLBI work. The Arecibo antenna has hardware features which make it equivalent to an 
azimuth-elevation mount. The Nancay array has been treated by Ortega-Molina (1985), but the 
model is not presently incorporated in MODEST code. 
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Figure 6. Schematic representations of the four major antenna geometries used in VLBI 


2.8.4 SITE VECTORS 

In the modeling software is the facility to provide a time-invariant offset vector in local 
geodetic coordinates (east, north, and local geodetic vertical) from this point (antenna location) to 
a point elsewhere, such as a benchmark on the ground. This is particularly useful in work involving 
transportable antennas which may be placed in slightly different places relative to an Earth-fixed 
benchmark each time a site is reoccupied. In modeling that offset vector, we make the assumption of 
a plane tangent to the geoid at the reference benchmark and assume that the local geodetic vertical 
for the antenna is parallel to that for the benchmark. With these assumptions there is an identity 
in the adjustments of antenna location with changes derived for the benchmark location. The error 
introduced by these assumptions in a baseline adjustment is approximately A B x (d/R^), where 
AS is the baseline adjustment from its a priori value, d is the separation of the antenna from the 
benchmark, and Re is the radius of the Earth. To keep this error smaller than 0.01 cm for baseline 
adjustments on the order of 1 meter, d < 600 meters is required. 
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More troublesome is that an error in obtaining the local vertical by an angle 5©, when using 
an antenna whose intersection of axes is a distance, H , above the ground, can cause an error of 
HtinSB « H6B in measuring the baseline to the benchmark (Allen, 1982). Unless this error is al- 
ready absorbed into the actual measurement of the offset vector, care must be taken in setting up 
the antenna so as to make SQ minimal. For a baseline error < 0.1 cm, and an antenna height of 
10 meters, SB < 20 arcseconds is required. Often plumb bobs are used to locate the antenna po- 
sition relative to a mark on the ground. This mark is, in turn, surveyed to the benchmark. Even 
the difference in geodetic vertical from the vertical defined by the plumb bob may be as large as 
1 arc minute, thus potentially causing an error of 0.3 cm for antennas of height 10 meters. Conse- 
quently, great care must be taken in these measurements, particularly if the site is to be repeatedly 
occupied by antennas of different sizes. 

2.8.5 FEED ROTATION 

Another physical effect related to antenna structures is the differential feed rotation for 
circularly polarized receivers, Liewer (1985) has calculated the phase shift 0 for various antenna 
types. It is zero for equatorially mounted antennas. For altazimuth mounts, 

tan 0 = cos <j> sin A/(sin 4> cos 6 — cos <f> sin 6 cos h) (2.211) 

with 4> = station latitude, h = hour angle, and S = declination of the source. For X-Y mounts, two 
cases are distinguished: orientation N - S or E - W. The respective rotation angles are 

tan(— 0) = sin 4> sin hj (cos <f> cos 6 ■+ sin <f> sin S cos h) {N S) (2.212) 

tan(— 6) = — cos h/(sin 6 sin h) (E — W) (2.213) 

The effect cancels for group delay data, but can be significant for phase delay and delay rate data. 
The effect on phase delay is 

T= (02-01 )// (2.214) 

where / is the observing frequency and 0, the phase rotation at station The feed rotation correction 
is an optional part of the MODEST model. 

Finally, another small correction which accounts for the effect of orientation of hour angle- 
declination (HA-Dec) and X-Y antennas on the tropospheric path delay was considered by Jacobs 
(1988). Details are given in the troposphere Section, 4.6. 

2.8.6 THERMAL EXPANSION 

By analogy with the model for atmospheric loading in Sec. 2.4.4, diurnal variations of the 
temperature cause vertical displacements of the antenna reference point. These can amount to several 
mm for ordinary day-to-night temperature variations for large antennas. If VLBI data acquired during 
a variety of weather conditions are to be processed simultaneously, it may be important to account 
for the vertical motion of the reference point. 

A rudimentary model of this effect is implemented in MODEST. It assumes that the vertical 
displacement A r of the antenna reference point, a distance h above the ground, is 


Ar = a(T — T re f)h (2.215) 

where ot is the coefficient of thermal expansion and T r «y is the reference temperature. The height 
h may be approximated as one half of the antenna diameter. The reference temperature is taken 
to be equal to the global average temperature at each station, or the universal average 292 K (used 
in tropospheric mapping in Sec. 4.2) if the former is not available. For steel and concrete, a linear 
expansion coefficient of 12 ppm is appropriate. Thus for a 70-m antenna, the vertical motion is 0.42 

mm/K. 
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2.8.7 ANTENNA SUBREFLECTOR FOCUSING 

For some experiments (notably the Time and Earth Motion Precision Observation [TEMPO] 
project) which use the DSN Cassegrain antennas, the subreflector is moved in order to maintain focus 
and to optimize the signal-to-noise ratio. Such motions introduce systematic errors in the antenna 
position derived from the measurements. For experiments performed in this “slew” mode, the path 
delay may change by «8 cm over the 6° - 90° elevation range. Simulations show that this effect is 
almost entirely absorbed by the clock epoch and local station vertical coordinate parameters. For 
baselines between two 70-m antennas, this causes a potential error of up to 12 cm in length. This 
effect is most easily modeled as a site vector relating fixed and slewed antenna positions. 
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2.9 PARTIAL DERIVATIVES OF DELAY WITH RESPECT TO 
GEOMETRIC MODEL PARAMETERS 

With respect to any given parameter, the calculation of the time-delay model must be at least 
as accurate as the data is sensitive to that parameter. Consequently, such effects as the curvature of 
the wave fronts were considered. However, such detail is not necessary for determining the derivatives 
with respect to the relevant model parameters. Here, the plane wave approximation is sufficient. 
Iteration on the estimated parameters and the rapid convergence of an expansion of the time delay in 
the relevant parameters about some a priori point permit this simplification. 

In this plane wave approximation we wish to obtain the parameter derivatives with respect to: 

1. the nominal baseline components (actually, station locations), 

2. the parameters of the whole Earth orientation matrix Q described in Section 2.6, 

3. the solid-Earth tidal parameters, 

4. the parameters of source location (right ascension and declination), 

5. the antenna axis offsets, 

6. the constant, in the retardation of the light ray due to gravitational effects. 

The expressions for these derivatives are considerably simplified if tensor notation, with the 
Einstein summation convention, is employed. Before proceeding any further, we make the following 
definitions for this section: 

t = time delay modeled in the geocentric frame, 

t • = this same time delay, but modeled in the Solar System center of mass frame, 
s = source unit vector (in the celestial system at rest with respect to the Solar System 
center of mass), 

P = velocity of the geocentric frame as measured in the Solar System center of mass frame 
(remember, all distances are measured in time; thus, this quantity is dimensionless), 

02 = velocity of station #2 in Solar System center of mass frame, 

p = 1 + 8 0 2 . This is a factor cs 1.0001, which arises from the motion of station #2 during 
the passage of the wave front from station #1 to #2, 

TT =(1-/?T 1/2 , 

•to 

Q = matrix which transforms from the terrestrial system to the celestial system, 

Lo = the baseline vector in the terrestrial system, 

L ( = this same baseline vector in the celestial system center of mass frame, 

L = this same baseline vector in the celestial system. 

With these definitions (2.177) may be written 



t = L. 

(2.216) 

For plane waves from (2.2): 

r _ & • [r 2 - ri] b L, b L, 

(2.217) 


1 — y 02 1 + 6 02 p 

Thus, 

'“-Til 

(2.218) 

For parameters (represented symbolically by r?) associated with L tk only: 



£— fa -**>7 + ■>*]*£ 

(2.219) 

Define the vector: 

+7 Pk] 

(2.220) 
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Then 


dr 

dr) 


Vk 


BL th 

dr) 


2.9.1 SOURCE PARAMETERS 

For parameters associated with the source position only: 


Since 


Define the vector: 
Then, 

For example: 
Then, 

and 

and 


I 1 = -7(1 -02. A) — 

Of) p 


d*k _ dp 
dr) p dr) 


P = 1 + 2 t 


dr] p 

--7(1 "A, ft) — 

P 


d*k _ Sk02i 
dr) p dr) 

3 k& 2 , 


Ski- 


ds) 

dr) 


Mj = -7(1 


Sii- 


&i&2, 


ds 

da 


_ as ds > 
dr) J dr) 

s = [ cos 6 cos a, cos 6 sin a, sin 6 } 

= [ — cos 6 sin a, cos 6 cos a, 0 J = [ A?, A3 ] 


ds 

— = [ — sin S cos a, — sin 6 sin a, cos 6 ] = [ Fi , F^ y F$ 
do 


dr 

da 

dr 

dS 


= MiAi 

= M X F X 


Or, if we define the matrices: 

and 

then: 


( Al M 

G = \ A 2 F 2 I 
v>h F 3 J 

M = ( M u M 2 ,M z ) 


dr dr 
da ’ d6 


= MG 


( 2 . 221 ) 


( 2 . 222 ) 

(2.223) 

(2.224) 

(2.225) 

(2.226) 

(2.227) 

(2.228) 

(2.229) 

(2.230) 

(2.231) 

(2.232) 

(2.233) 

(2.234) 
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For a linear model of source “proper motion" [Eqs. (2.9l)-(2.92)], the partials of r with respect to 
the time rates of change of right ascension and declination (d, 6) are 


dr dr 
da 1 36 


(t-t 0 )MG 


(2.235) 


where to is a reference time. 

The partials of the observables with respect to the flux ratio K and the baseline projection R for 
a two-component model are obtained from Chariot’s (1990a) expression for the group delay: 

dr _ 2nR\\ - cos(2tt.R)] K* - 2 K* - 2 JC 2 -2K+1- 4JT 2 cos(2 kR) 
dK W (1 + Kf [jp + 2Xcos(2irJl) + l] : 2 

dr _ 2irK(l - K ) 

dR~ w(l + K) 

(1 + + 2xR ( 1 + K) sin(2ff.R) — [(l — jjf) 2 -f 2/T cos(2ri?)| cos(2?rj?) 

[K ‘ 2 + 2K cos(2 nR) + l] 2 

Finally, the partials with respect to the remaining two structure parameters s and 6 are obtained from 
dr/dR via 


(2.236) 


(2.237) 


dr dr dR . _ dr 

Ts = dRj7= (“« n * + «cos 6) — (2.238) 

dr dr dR dr 

de = dRle =a{ucos6 ~ vsine) dR ( 2 - 239 ) 


2.9.2 STATION PARAMETERS 

For station location parameters the algebra is somewhat more complex. Since 

L * -r 2 {h) -ri(ti) (2.240) 

= r 2 (t 2 ) -ri(ti) - 0 3 [t 2 ~ty] 

= r 2 (t 2 ) -r x (t!) -70 2 !r 2 (ti) -r'(ti)] 0 

= + h ~ l)[»i(*i) - ri - ifi 2 fi |r' (*') - ,'(*')] 


we have: 

or in tensor notation 

Define the tensor: 

Then 

Since 


L, = L + (7 - 1)L fifi - if} 2 /3 L 



(1 - m 

d 7 


- 102. 



E<j = <5tj + 


(l ~ l)ft 

0 2 


- 102 , 


0 : 


L tl — EijLj 

Lj = Q : kLoi, 

L., = E XJ QjkL 0k 


(2.241) 

(2.242) 

(2.243) 

(2.244) 

(2.245) 

(2.246) 
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T1 "“’ r - .£o. 

For par m« which are involved with station locations expressed in thr terrestrial coordin.tr system: 




(2.248) 


where the vector element (2.249) 

Bk = 'Vit’ijWjk 

Such parameters are: rj,., (radius oK spin axis), A? (longitude), a? (height above the equator) t 

refers to station number, i.e., i = 1, 2. Define the matrix. 

W = [-R U R 3 , -Ax, A 2i -Zi, Z 2 , -ft , * 2 , -Ai, A 2i -Zi, Z 2 , -ft, V 2| -ft, ft, -*i, *a] ( 2 25 °) 

where each column contains the partials of the I 0 component vectors x, y, z with re^ct to the 
parameters. For example, for the constant terms in the cylindrical station coordinates [see Eq . 

(2.36) through (2.38)]: 


ft = 


f dLo, ^ 

Wp, 

dLo, 



A, = 



dr°.r. 

dLo. 

\ Wp, J 

( dLo. \ 

SA? 
dL o v 
SA° 

dLo. 

V 3a° ; 

( dLo± \ 

dz° 

dLo, 

dz? 

dL 0 . 

\WJ 

For the station coordinate rates, 

ft = (t - t 0 )ft A, = (t-t 0 )A, Z, = (t - t 0 )Z, 
From Eqs. (2.48) through (2.59), and relying on Williams (1970): 


(2.251) 


(2.252) 


Z. = 


(2.253) 


(2.254) 


V, = 


( \ 

dh 2 i 

dt oy 

dh 2 , 

db xz 

V ah 2 , ) 


= S(i)V(i)W(i 



(2.255) 
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Hi = 


( dkz \ 

dhi 

d6iy 

~dhl 

dS xx 

\ dl 2 < J 

( dSix \ 

dipi 


0 

,< 2 ) ( 


= 5(t)V(,)W(t) g\ (0 

U 3 , M. 


(2.256) 


4>. = 


*v 


£4 

dipi 

dS xx 

V dip { J 


= S(i)V(i)W(i) 


( gg ( i 2) (0 A 

di>i 

dg^QQ 

9A 

dg 3 2> (0 

V 5l6, J 


(2.257) 


where » = 1 implies station #1, i = 2 implies #2, and S(l) = —1, while S(2) = 1. These partials of 
with respect to ip are 


= ^lhr p • K.\y P X. - x p Y,} 


dip 

a„(2) 




dip 

Rl 

* J2) 

d 93i _ 

_ 3m. r\ 


dip 

Also, define a vector: 


yj*l + yp l 

£ 5 P ^ [VP^* “ X P^* 


jr p • R.HipX, + y p y.] - [ipY. - y P X.] 2 


y/W*- / = l 2r p R * - *p-Z,] 

V x p ^ 2 


(2.258) 

(2.259) 

(2.260) 


D = 


dr 

dr dr dr dr 

dr 

dr 

dr 

dr 

dr 

dr 

dr 


’ dr 0 .p, ’ SA°’ 5A°’ dz\ 

' dz°' 

dr 4px 

’ ar. n ’ 

dXi 1 

5A 2 ’ 

dz\ 

di 

dr 

dr dr dr dr 

dr 

dr 

dr 






dTi 1 

562 1 dhTi’ 5^7’ aiZ' 

1 w 
« 
— j 

|<0 

Wi' 

dip 2 







Then 


D = BW 


(2.261) 

(2.262) 


2. 9. 2.1 Ocean Loading Parameters 

Partials with respect to the amplitudes [Eq. (2.74)] of the ocean loading model are 
trivial in the local coordinate system. Transformation to geocentric coordinates via (2.68) yields the 

dr 

corresponding partials with respect to the amplitudes, . 

2. 9. 2. 2 Atmosphere Loading Parameter 

Partials of the station vertical displacement with respect to the atmosphere loading scale 

factor /, = — 0.55£, are transformed via (2.68) to give the observable partials — j. 

d f of 
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2.9.3 EARTH ORIENTATION PARAMETERS 

Certain parameters such as UTl, polar motion, precession, and nutation affect Q only. For 
these parameters, symbolized by tj, 


Define a vector: 




Ki = * k E ki 


§r* ^ >»• 


(2.263) 


(2.264) 


(2.265) 


for parameters which affect only the orientation of the Earth as a whole. 

2.9.3. 1 UT1 and Polar Motion 

A number of parameter partials are available for the orientation of the Earth. These are for 
UTl, X pole, Y pole, and nutation, as well as the angular offset and angular rate terms in the Earth 
orientation perturbation matrix fi. From (2.98): 


Define the matrix: 


Q = n PNUXY 


SY I 

Y' = = I 0 — sin ©2 cos ©2 

^©2 l o -cos ©2 —sin ©2 


Then, the partial required for the Y polar motion parameter is: 


22 - = n PNUXY' 
d© 2 


(2.266) 


(2.267) 


(2.268) 


dX 

An analogous technique is used for the X pole angle, working with the matrix partials Partials 

with respect to UTl involve a slight complication due to the last three terms in Eq. (2.105). On the 
assumption that only the term linear in T u contributes significantly, 

— - = |^.(l + l/365.25) (2.269) 

d(UTl) dH x 


Partials of time rates of the UTPM parameters are related to those of the parameters themselves: 


dQ dQ U t\ 


(2.270) 


2.9.3. 1.1 UTPM Tidal Amplitudes 

Partials of r with respect to the amplitudes of both the solid and ocean tide contributions 
to UTl and polar motion are also easily obtained. Again, the required partials of Q are written as 

^ = nPN^-XY ( 2 . 271 ) 
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and 


^ = CIPNU^-Y 

or) or} 


with a similar expression for the Y component. For the solid tide amplitudes A< (2.113), 


dU 


du ae 3 


dA, 9{UT1) dAi 


(2.272) 


(2.273) 


and 


9©3 

JZ 


= tin 


L >=i 


For the ocean tidal amplitudes A f and B? of (2.114), 


ax? 


= cos 


^ j a j + ^)| 

V=i J 


(2.274) 


(2.275) 


2.0. 3.2 Nutation 

Partial derivatives of the VLBI observables with respect to the nutation angles and amplitudes 
appear formidable at first sight, but are relatively easy to evaluate if the calculation is performed in 
an organized fashion. Symbolizing the parameters by 17 , we need to evaluate the partials of the matrix 
Q with respect to rj: 

^ = np (^ u+N i£L)^r LxY (2 - 276) 

dr) \96\p 96ip J dr) 


dr) \ d6e 96c J drj 


(2.277) 


dN 


Since 6 c = e-£, the first partial on the right hand side of Eq. (2.277) is equal to The derivatives 
of N with respect to the angles 6ifi and 6c are easily obtained from the expression for N in Eq. ( 2 . 120 ). 


dN 


— sin Sip 


: cos 6rj) 


sin £ cos 6\l> 


= — cos £ cos 6ip — cos £ cos 5 sin 6ip — cos £sin esin Sip 

\ _ sin £ cos 6xp - sin £ cose sin 6$ - sin £sin esin 6\l> 


(2.278) 


and 


dN 


— sin e sin 


cos c sin 6rp 


1 - 0 - cos £ sin e cos <5 \p + sin £ cose cos £ cos £ cos Sip + sin i sins 

' 0 - sin £ sine cos 6$ — cos £ cos £ sin£ cos £ cos Sip - cos £sin e 


From Eq. (2.103), the partials of U with respect to 6\p and 6e are 


dU 


d5ip } Sc 


- sin H -co sH 

cos H — sin H 
0 0 


dH 

96\p % 6c 


(2.279) 


( 2 . 280 ) 


and, from Eq. (2.108), 


= cose / (cos 2 6ip -r cos 2 c sin 2 <5t/>) 

96 xP 1 V 

— = - sine tan 6$ / (1 + cos 2 e tan 2 Sip) 
d5e 


( 2 . 281 ) 

( 2 . 282 ) 
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the [^-dependent terms in Eqs. (2.276) and (2.277) are evaluated. 

Partials of Sip and Se with respect to the parameters A x} and Bn are obtained immediately from 

Eqs. (2.130-2.133). For the “free nutations”, 




dSiP J . _ 

= sm w/T, 

OAqo 

dSt 1 
dB oo 

cos a;/ T 

(2.283) 

d6rf> J rp . T 

-Tsmw/T, 

oAio 

d6e f __ 
35io 

T cosw/T 

(2.284) 

d6ij> f T 

dA 20 

36c' 
dB? o 

sin 

(2.285) 

_ T 
= T cos 1 , 

dA 30 

36c' _ 
dB 30 

T sin w/T 

(2.286) 

and for the 1980 IAU series terms (j = 1 to 106): 





36e 
3 Boy 

5 

cos [5Z 
»=1 

(2.287) 

dA 1} . 

3<5e 

3B iy “ 

5 

= T cos [T: fcyiai(T) 1 

«=i 

| (2.288) 


d6e 

dB 2] " 

5 

»=1 

(2.289) 

|^ = Tcos[£fey,a,(T)] , 

d6e 

BB 3 y ‘ 

5 

= Tsin[^^ai(T) 

»=i 

] . (2.290) 


2.9.S.S Precession 


Partial derivatives of the observables with respect to precession parameters are evaluated 
in a manner similar to those with respect to nutations. Symbolizing either the luni-solar precession 
Pls or the planetary precession ppl by the partial of the rotation matrix Q is 



(2.291) 


The partials 
element Pa- 


dp 

d7T 


are very complicated, and will be written in terms of the partials of each matrix 


— = — cos £o sin £ cos © cos Z / 2 — sin £o cos C sin ® cos ^ 

T dpLS 

— cos Iq cos f cos © sin Z/2 — cos £q cos £ sin Z/2 

— cos £o sin $ cos Z/2 (2.292) 


— — sin c cos © cos Z/2 H- cos £ cos © sin Z/2 

T dppi 

+ cos f sin Z/2 -r sin <; cos Z/2 (2.293) 


1 dP 

— = - cos Co sin f cos © sin Z/2 - sin e 0 cos f sin © sin Z 

T dpLs 

+ cos €q cos £ cos © cos Z/2 -h cos fo COS f cos Z/2 

- cos sin r sin Z/2 (2.294) 
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(2.295) 


1 d P 

— *r — — = sin $ cos © sin Z/2 — cos $ cos © cos Z/2 
T OPPL 

— cos f cos Z/2 4 sin £ sin Z/2 


j_£Pl3 

r dpLs 


— cos £o sin $ sin © / 2 4 sin £q cos £ cos © 


^ — = T sin £ sin ©/2 


1 dP 2l 
T dpLS 


— cos io cos £ cos © cos Z/2 4 sin £o sin f sin © cos Z 

4 cos io sin £ cos © sin Z/2 4 cos £o sin £ sin Z/2 

— cos £q cos £ cos Z/2 


1 3P 2l 
T SprL 


= cos £ cos © cos Z/2 — sin £ cos © sin Z/2 
- sin £ sin Z/2 4 cos ; cos Z/2 


1 dP22 
T dpLS 


cos So cos f cos © sin Z/2 + sin £o sin £ sin © sin Z 

cos £ 0 sin £ cos © cos Z/2 — cos £o sin £ cos Z/2 
cos £q cos £ sin Z/2 


1 <9P 2 2 

— t = cos f cos © sin Z/2 4 sin £ cos © cos Z/2 

T opp L 

4 sin £ cos Z/2 4 cos £ sin Z/2 


1 5/^3 
3T 3pz,s 


1 dPzi 
r dpz,s 


1 5^32 
T 5p£,c 


— cos £o cos £ sin 0/2 — sin £o sin £ cos © 

- -- 3 = T cos £ sin ©/2 

dpri, 

— sin £o cos © cos Z 4 cos £o sin © sin Z/2 

= -T sinG sin Z/2 

dpPL 

— sin ?o cos © sin i? — cos £q sin © cos Z/2 


dP- 


Si 


dppL 

dP 33 

dpLS 


= T sin © cos Z/2 
= — T sin fo sin © 


dP- 


33 


dp pl 


= 0 


A check on the algebra may be performed by noting that 


d_P_ 

dppL 




(2.296) 


(2.297) 


(2.298) 


(2.299) 

(2.300) 

(2.301) 

(2.302) 

(2.303) 

(2.304) 

(2.305) 

(2.306) 

(2.307) 

(2.308) 

(2.309) 


54 



and 


dP 

dpLS 


dP _ . dP 

= - cos l 0 - + Tsinfo^s 

dp pl oW 


The corresponding partials of the V matrix are much simpler: 

( — sin H — cos H 

= T cos £q 


dU 


dpLS 



au_ 

dppL 




(2.310) 

(2.311) 

(2.312) 


2.9. 8.4 Rotational Tweaks 


Finally, the partials of the nutation matrix with respect to the “tweaks" Ai/i and Ac are 

dN dN 

obtained by making the replacements (2.134) and (2.135) in N. The partials an ° are 

then seen to be identical to Eqs. (2.278) and (2.279), with the same replacements for Sip and Sc. 
Expressions analogous to Eqs. (2.280-2.282) account for the shift of the equinox by nutation changes 
Sip and Sc. If the a priori tweaks are zero, the partials are exactly equal to the expressions (2.278) 
and (2.279). 

For the parameters in the perturbation matrix, ft, from (2.165-2.166): 


(2.313) 



(2.314) 


where £ is the number of years from the reference epoch [e.g. t J2000). Then, by substituting these 
matrices for ft in (2.98), we obtain the appropriate partials of Q for perturbations about the x axis. 
By analogy, the perturbation parameters about the y and z axes may also readily be obtained. 

2.9.4 ADDITIVE PARAMETERS 


If we seek the partials of parameters that affect only the “add-on" terms in r = tq + At, then 
from (2.177) we have: 


drj 

for terms which were “added on" in the Solar System barycenter. An example is gravitational bending: 

dr „ v Ac 


For terms “added on” in the geocentric frame, then: 

dr _ dA r 
drj dr] 

An example is the antenna offset vector. In this case: 


(I + 'Ww) 


(2.315) 
ending: 

(2.316) 


dr 


and 


<9(offset station #2) 
dr 


±\A-P- 


i ] 2 


= ±J\ - (s I' 2 


(2.317) 


(2.318) 


(2.319) 


^(offset station #1) 

where the choice of sign for each station is determined by the choice of sign for that station in the 
model portion. 


55 



SECTION S 


CLOCK MODEL 


The frequency standards (“clocks”) at each of the two antennas are normally independent of each 
other. Attempts are made to synchronize them before an experiment by conventional synchronization 
techniques, but these techniques are accurate to only a few /xs in epoch and ft 10“ 12 in rate. More 
importantly, clocks often exhibit “jumps” and instabilities at a level that would greatly degrade 
interferometer accuracy. To account for these clock effects, an additional “delay” r c is included in the 
model delay, a delay that models the behavior of a station clock as a piecewise quadratic function of 
time throughout an observing session. Usually, however, we use only the linear portion of this model. 
For each station this clock model is given by: 


T c — Tel + T C 2[t tref] + ^re/) /2 


(3.1) 


The term, t re /, may be set by the user as a specific time (Julian date), or by default taken as the 
midpoint of the interval over which the a priori clock parameters, r cl , r c 2, r c 3, apply. 

In addition to the effects of the lack of synchronization of clocks between stations, there are 
other differential instrumental effects which may contribute to the observed delay. In general, it is 
adequate to model these effects as if they were “clocklike”. However, the instrumental effects on delays 
measured using the multifrequency bandwidth synthesis technique (Thomas, 1981) may be different 
from the instrumental effects on delays obtained from phase measurements at a single frequency. 

The bandwidth synthesis process obtains group delay from the slope of phase versus frequency 

(r = — ) across multiple frequency segments spanning the receiver passband. Thus, any frequency- 
independent instrumental contribution to the measured interferometer phase has no effect on the 
delay determined by the bandwidth synthesis technique. However, if delay is obtained directly from 

the phase measurement, ^>, at a given frequency, v, then this derived phase delay (r p d = — ^ 
have that instrumental contribution. 

Because of this difference, it is necessary to augment the “clock” model for phase delay measure- 
ments: 

— r c + r c 4(t tre/) "b T e &{t £ re /) 2 /2 (^*2) 


does 


where r c is the clock model for bandwidth synthesis observations and is defined in (3.1). Since 
the present system measures both bandwidth synthesis delay and phase delay rate, all of the clock 
parameters described above must be used. However, in a “perfectly” calibrated interferometer, 
r c 4 = r c s = 0. This particular model implementation allows simultaneous use of delay rate data 
derived from phase delay, with delay data derived by means of the bandwidth synthesis technique. 
However, our particular software implementation currently is inconsistent with the simultaneous use 
of delay derived from bandwidth synthesis and delay obtained from phase delay measurements. 

An optional refinement of the clock model is also available. For dual-frequency (S/X) delays, an 
additional clock offset may be estimated. It originates from the differential instrumental delay and 
fringe fitting delays for S- and X-band data, which may be sizeable. For dual-frequency observables, 
the clock model depends on this differential instrumental delay and on the frequencies a>5,u>x in the 
individual bands as 

r c6 ^x/t^x “*^1) (3.3) 


The differential instrumental delay r cG is normally highly correlated with the usual clock offset r cl , 
but under some circumstances may convey additional information. 

To model the interferometer delay on a given baseline, a difference of station clock terms is 
formed: 

elation. 2 T C, ta ,t lon j (3.4) 

Specification of a reference clock is unnecessary until the parameter adjustment step, and need not 
concern us in the description of the model. 

The partial derivatives of model delay with respect to the set of six parameters (r cl , r c2 , r c3 , r c4l t c5 , 
and r c c) for each station are so trivial as to need no further explanation. 
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SECTION 4 


TROPOSPHERE MODEL 

In order to reach each antenna, the radio wave front must pass through the Earth’s atmosphere. 
This atmosphere is made up of two components: the neutral atmosphere and the ionosphere. In 
turn, the neutral atmosphere is composed of two major constituents: the dry and the wet. The dry 
portion, primarily oxygen and nitrogen, is very nearly in hydrostatic equilibrium, and its effects can be 
accurately estimated simply by measuring the barometric pressure. Typically, at sea level in the local 
zenith direction, the additional delay that the incoming signal experiences due to the troposphere is 
approximately 2 meters. Except for winds aloft, unusually strong lee waves behind mountains (e.g., 
Owens Valley, California), or very high pressure gradients, an azimuthally symmetric model based on 
measurements of surface barometric pressure is considered adequate. We have not yet investigated 
where this assumption breaks down, though “back-of-the-envelope" calculations indicate that, except 
in the unusual cases above, the error in such an assumption causes a less than 1-cm error in the 
baseline. 

Unfortunately, the wet component of the atmosphere (both water vapor and condensed water in 
the form of clouds) is not so easily modeled. The experimental evidence (Resch, 1983) is that it is 
“clumpy", and not azimuthally symmetric about the local vertical at a level which can cause many 
centimeters of error in a baseline measurement. Furthermore, because of incomplete mixing, surface 
measurements are inadequate in estimating this contribution which even at zenith can reach 20 to 
30 cm. Ideally, this tropospheric induced delay should be determined experimentally at each site. 
This is particularly true for short and intermediate (B < 1000 km) baselines, where the elevation 
angles of the two antennas are highly correlated in the observations. For long baselines, both the 
independence of the elevation angles at the two antennas, and the fact that often the mutual visibility 
requirements of VLBI constrain the antennas to look only in certain azimuthal sectors, allow the use 
of the interferometer data itself to estimate the effect of the water vapor as part of the parameter 
estimation process. For this reason, and because state-of-the-art water-vapor measurements are not 
always available, we also have the capability to model the neutral atmosphere at each station as a 
two-component effect, with each component being an azimuthally symmetric function of the local 
geodetic elevation angle. 

At each station the delay experienced by the incoming signal due to the troposphere can be 
modeled using a spherical-shell troposphere consisting of a wet component and a dry component: 

Ttrop station t = T%vet irop ^dry trop (4.l) 


The total troposphere model for a given baseline is then: 


Tt Ttrop station 2 T trop station 1 (4.2) 

If Ei is the apparent geodetic elevation angle of the observed source at station i, we have (dropping 
the subscript »): 

T trop = PZ+ry Rdry [&) + PZ wt * ^wet (^) (4.3) 

where pz is the additional delay at local zenith due to the presence of the troposphere, and R is an 
elevation angle mapping function. 

For some geodetic experiments, the observed delay has been corrected for the total tropospheric 
delays at the two stations, which are in turn calculated on the basis of surface pressure measurements 
for the dry component, and water-vapor radiometer measurements for the wet component. This 
correction is recorded in the input data stream in such a way that it can be replaced by a new model. 
In the absence of such external calibrations, it was found that modeling the zenith delay as a linear 
function of time improves troposphere modeling considerably. The dry and wet zenitli parameters are 
written as 

pz*.„ = P° Zi ,„ + PZ d .„{t- to) (4.4) 
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where to is a reference time. 

Since the model is linear in the parameters p° and p, the partial derivatives with respect to zenith 
delays and rates are trivial. They are: 


dr 




— «r w 


(4.5) 


and 


dr 


w 


(t ^o)/(0*®»< PT w 


where /(t) = 1 for station #2, and —1 for station #1. 


4.1 CHAO MAPPING FUNCTION 


(4.6) 


The simplest mapping function implemented in MODEST code is that obtained by C. C. Chao 
(1974) through analytic fits to ray tracing, a function which he claims is accurate to the level of 1% 
at 6° elevation angle and becomes much more accurate at higher elevation angles. 


R = 


1 


sinE + 


A 

tan E -b B 


(4.7) 


where 


A dry = 0.00143 (4.8) 

B drv = 0.0445 (4.9) 

A wet = 0.00035 (4.10) 

B wet = 0.017 (4.11) 


The user must specify values for the zenith delays. 

The partial derivatives of delay with respect to the parameters A dry and B dry are: 

3 ^ = -f(i)pz^Rl v /(t^E + B drv ) 


(4.12) 


and 


dr 

dBd ry 


f[})pZdr y ^dry-Adry / E ■+■ Bdry) 


where Rdry is the Chao mapping function, and E is the elevation angle. 


4.2 LANYI MAPPING FUNCTION 


(4.13) 


Analyses of intercontinental data indicate that the Chao mapping function [Eq. (4.7)] is inade- 
quate. To rectify this situation, two modifications have been made to the MODEST code. First, the 
dry-troposphere mapping parameters Adry and Bdry of the Chao mapping function Rdry have been 
promoted to the status of estimable parameters. Second, the code now permits the use of two more 
accurate mapping functions. The first of these is the analytic function developed by Lanyi (1984). In 
its simplest form, this mapping function employs average values of atmospheric constants. Provision is 
made for specifying surface meteorological data acquired at the time of the VLBI experiments, which 
may override the average values. Using numerical fits to ray-tracing results, Davis ct al. (1985) have 
arrived at another function, designated the CfA-2.2 mapping function. Comparisons indicate that the 
Lanyi and CfA functions are in agreement to better than 1 cm over an extreme range of atmospheric 
conditions down to 6° elevation angles. Finally, an approximate partial derivative is obtained with 
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respect to one parameter in the Lanyi mapping function; this permits adjustment even m the absence 
of surface data. The Lanyi function was made the default MODEST troposphere model in ear y 1986. 

Motivation for and full details of the development of a new tropospheric mapping function are 
given by Lanyi (1984). Here we attempt to give a minimal summary of the final formulas. I he 
tropospheric delay is written as: - 

itrop = F{E)/ sin E (4-14) 

where 

F{E) = Pz dry F drv (E) + Pz m .,F W et{E) 

+ I Pz*. y Fbi(E) + 2 P z irt Pz„.,F M (E) + p 7 z m . t Fba[E)}/& + p\ d . y Fb*[E) / A (4.15) 

The quantities PZdry and pz^ have the usual meaning: senith dry and wet tropospheric delays. 
A is the atmospheric scale height, A = kT 0 /mg c , k = Boltzmann’s constant T 0 = average surface 
temperature, m = mean molecular mass of dry air, and = gravitational acceleration at the center of 
gravity of the air column. With the standard values k = 1.38066 x 10 «g/K, "* “ 4 8097 X 10 . *' 

9c = 978.37 cm/s 2 , and the average temperature for DSN stations T 0 = 292 K, the scale height 

The dry, wet, and bending contributions to the delay, F dr y[E), ftoet(-^)i an< ^ f6i,b2,63,M(-E)> are 
expressed in terms of moments of the refractivity as 

Fdr V (E) = A lo (E)G(AAfuo, u) + 3£ruM 2 ioG 3 (Mno. «)/2 
F wet (E) = A 0 i(£)G(AMioi/Afooi,u)/M>oi 
Fbi{E) = [(7G 3 (Mno>u)/sin 2 iP- A/ o 2 oG 3 (Mi 2 o/Md 20. «)]/2tan 2 £ 

F b2 (E) = -AfoiiG 3 (M lu /Moii,u)/2Mooitan 2 E 

Fbz[E) = — Mo02G 3 (Mi02/A^002i u )/2Mo 01 tan2 E 

Fb4{E) = M 0 3oG 3 (M 13 o/-^03o,u)/tan 4 E 

A misprinted sign in the last of Eqs. (5) of Appendix B of Lanyi (1984) has been corrected in 
Eq. (4.21). Here G(g,u) is a geometric factor given by 

G[q, u) = (1 + gu) -1 ^ 2 (4-22) 

WltH u = 2cx/ tan 2 E (4-23) 

where a = A/i2 is a measure of the curvature of the Earth’s surface with standard value 0.001345. 

The quantities A, m (£) and M, im are related to moments of the atmospheric refractivity, and are 
defined below. A 10 (E) involves the dry refractivity, while A 0 i{E) is the corresponding wet quantity. 
The Aim(E) are given by 


(4.16) 

(4.17) 
(4-18) 
(4.19) 
(4-20) 
(4-21) 


^ " ( — l)"+*(2n— l)!!M n - fc ,i.r 

Aim[E) — Maim + / , 2 2 rl Jt'(n — Ik)’ 

n=lk = 0 ‘ V 


■k,l,m j 

ti 

n 

A Milm 


1 -f AuAfiim/Moim 


A/q Im 


(4.24) 


with the scale factor A = 3 for E < 10° and A = 1 for E > 10°. Only the two combinations 
(l, m) = (0,1) and (1,0) are needed for the A im {E). The moments of the dry and wet refractivities 

are defined as M 

(4.25) 


M n „ = J dq q n f l drv (q)fl,ct{<}) 


where f drv „, et (<j) are the surface-normalized refractivities. Here, n ranges from 0 to 10, i from 0 to 
3, and j from 0 to 2; not all combinations are needed. Carrying out the integration in Eq. (4.25) for 
a' three-section temperature profile gives an expression for the general moment M ni , : 


M n „ In' = {l-r a ")/a 


n+l + e -o<l 


‘[i-T'- +n+i ( 9 i, 9 2)]n 

i = 0 


b + i + 1 


& C 


- a2 ‘^ + " + i(gi,92)/a" + i 


( 4 . 26 ) 
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Here, 


72(<7i> Qt) = 1 — [92 - 


(4.27) 


The quantities qi and q 2 are the scale-height normalized inversion and tropopause altitudes, respec- 
tively. For the standard atmospheric model, qi = 0.1459 and q 2 = 1.424. The constants a and b are 
functions of the dry (a = 5.0) and wet (0 = 3.5) model parameters, as well as of the powers of the 
refractivities (i and j) in the moment definitions. Table IX gives the necessary a’s and b'&. 


Table EX 

Dependence of the Constants a and b 
on Tropospheric Model Parameters 


t 

; 

a 

b 

1 

0 

1 

a — 1 

0 

1 

p 

af) — 2 

2 

0 

2 

2(a-l) 

1 

1 

0+1 

0(a + 1) - 3 

0 

2 

u 

2(ap - 2) 

3 

0 

3 

3(a-l) 


Note that the normalisation is such that A /010 *= 1; this moment has therefore not been explicitly 
written in Eqs. (4.16) through (4.21). 

At present, provision is made for input of four meteorological parameters to override the default 
(average) values of the Lanyi model. These are: 1) the surface temperature To, which determines 
the atmosphere scale height; 2) the temperature lapse rate W t which determines the dry model 
parameter a; 3) the inversion altitude hj, which determines qx\ and 4) the tropopause altitude h 2 , 
which determines q 2 - A fifth parameter, the surface pressure po, is not used at present. Table 
X summarizes the four parameters and derived quantities, and their default values. Approximate 
sensitivities of the tropospheric delay (at 6° elevation) to the meteorological parameters are given in 
the last column. 


Table X 

Surface Meteorological Parameters in the Lanyi Mapping Function 


Parameter 

Default 

Derived parameter and value 

Sensitivity (6°) 

To 

292 K 

A = kTo/mg c = 8.567 km 

— 7 mm/K 

W 

6.8165 K/km 

q = 100To/WA = 5.0 

20 mm/K/km 

hi 

1.25 km 

9! = hi/ A = 0.1459 

—20 mm/km 

h 7 

12.2 km 

52 = h-2 j A = 1.424 

5 mm/km 
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Partials of the delay with respect to the dry and wet zenith delays are obtained from Eqs. (4.14) 
and (4.15): 


dpz*,, 


= /(0[W£) + 2p Zirv Fbi(E)/^]/smE 

+ [2 P z„ME)/A + *p\ lrv FMI* 2 ]l™E 

■ = m[F wei (E) + 2p Zi , t Fb2(E)IA + 2pz^F bi (E)/A}/smE 


t analvfiK Of data for which meteorological parameters are not available, it is convenient to introduce 
^approximation into the mapping function [Eqs. (4.14) and (4.15)1 which involves a on-P-ameter 
estimate. This parameter p accounts for deviations from standard meteorological conditions. T 

tropospheric delay is expressed as 


Ttrov = (PZ<rv+PZ—)/ BijlE + 1P 


dlfrop 


(4.30) 


where the partial derivative is 


dr trop _ [PZjn, 

dp ~ G(A/no,«)[l + G(Muo,u)] 8 inf; 

pz w .,u [Mll° ~ (4.31) 

+ G(Muo,u)G(M 1 oi/Mooi,u)[G(Wno,«) + G(M 1 oi/Mooi,«)]sinE 


4.3 CfA MAPPING FUNCTION 

Another approach to improved modeling of tropospheric delay was published by Davis et al. 
(1985). Analytic fits to ray-tracing results yield the CfA-2.2 mapping function 

p _ l (4-32) 


sin£ + £ 

tan E + — — — 
sin E + c 

where E is the elevation angle. The three parameters a, 6 , c are expressed in terms of meteorological 
data as 


a = 0.0002723 ( 1 + 2.642 x 10 _ 4 p 0 - 6.400 x 10 4 e 0 + 1-337 x 10 3 T 0 

- 8.550 x 10 _2 o - 2.456 x 10 _2 h 2 ) 

b = 0.0004703 ( 1 ■+• 2.832 x 10 -5 po + 6.799 x 10 -4 eo + 7.563 x 10 3 T 0 

- 7.390 x 10 -2 Of - 2.961 x 10 -2 h 2 ) 
c = - 0.0090 


(4.33) 

(4.34) 

(4.35) 


Here p 0 is the surface pressure and e 0 the surface partial water vapor pressure, both measured in 
millibars The quantities T 0 , a, and h 2 have the same meaning and units as in Section 4^2. Thu 
Tu„“r„ is one of five optional mapping functions in the MODEST model In connect, on wrth testing 
parameter estimation for the Lanyi function, the partial derivative of delay with respect to surface 
temperature T 0 in the CfA-2.2 function was also evaluated. It is 

d T PZdcyKjryl 3 - 641 x 10~ C ( sin£;+ c)|tanE-h fc/(sinE + c)] - 3.557 x 10 c a] ^ ^ 
37 ”” (sin E ■+■ c) [tan E -r 6 /(sin E 4- c)] 
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4.4 HERRING 1992 (MTT) MAPPING FUNCTION 

A further refinement of the CfA approach by Herring (1992) has yielded the following mapping 
functions for the dry and wet tropospheric delays, respectively: 


1 Q dry,wet 


Rdry ,wet — 


1 + 




lry,wct 


1 Cdry,wet 


sin E + 


&dry,wet 


b; 




(4.37) 


The parameters (a, fc, c)dry,wet depend on the site latitude (#), altitude (/f, km), and surface temper- 
ature (To, K) via 

a dry = 0.0006232 ( 1 + 0.02230 cos <f> - 0.03354 H + 0.003450T o ) (4.38) 

b d ry = 0.0025779 ( 1 - 0.06207 cos <f> - 0.01284# + 0.000799T o ) (4.39) 

Cdry = 0.071839 ( 1 - 0.05976 cos 4> - 0.00207# - 0.000029T 0 ) (4.40) 


a we i = 0.000187 ( 1 - 0.0590 cos <f> - 0.2787# + 0.00750T o ) (4.41) 

b wet = 0.000836 ( 1 -0.1221 cos <£- 0.1209# -f 0.00239T 0 ) (4.42) 

c wet = 0.04160 ( 1 - 0.0459 cos <f> - 0.0310 H + 0.00036T 0 ) (4.43) 

This function is becoming known in the literature as the “MTT" function. 

4.5 NIELL MAPPING FUNCTION 

Two new mapping functions (wet and dry) which use no surface meteorological data have 
recently been introduced by Niell (1994). They take the same continued fraction form as the Herring 
functions above (4.37). The parameters (a, b } c)dry,wct have a more complicated latitude and seasonal 
dependence, and there is also a (dry only) altitude correction. For latitudes within 15° of the poles, 
the parameters are constant, while in mid-latitudes they vary both with latitude and day of year (with 
proper adjustments for the reversal of seasons in the Southern Hemisphere). The altitude correction 
is proportional to altitude, and has an elevation dependence similar to that of (4.37). Because of the 
complicated functional dependences and numerous coefficients, the Niell mapping functions are not 
presented here. They are, however, one of the options in MODEST. 

4.6 ANTENNA AXIS OFFSET ALTITUDE CORRECTION 

Antennas with non-zero axis offsets, whose second rotation axis (A in Figure 5) moves vertically 
with changing orientation, have zenith troposphere delays that may vary by 1 to 2 mm. Equatorial and 
X-Y mounts fall in this class (see Figure 6). At low elevation angles this zenith variation is magnified by 
the mapping function to 1-2 cm. These variations must be modeled in experiments whose accuracies 
are at the millimeter level ( c.g . short-baseline phase delay measurements). Memoranda by Jacobs 
(1988, 1991) derive the corrections based on considering only the dry troposphere component, and 
include all terms necessary to achieve an accuracy of a few millimeters. The correction to be added 
to the zenith dry tropospheric delay is 


Sr = -pz 4 rv {H/A) \p (4.44) 

where H is the antenna axis offset, A the dry troposphere scale height (« 8.6 km), and ip is an angular 
factor that varies with the type of mount. For equatorial mounts, 

xp = cos <p cos/i (4.45) 
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where <t> is the geodetic latitude and h the local hour angle east of the meridian. The Richmond 
antenna correction has this form with <f> replaced by <$> W and h by a pseudo-hour angle h R (see Section 
2.8.3), where 

h K = arctan [cos Esin(0-«)/ [cos^iy sin E- sin cos Ecos(0 + 4 )l] t 4 ' 46 ) 

For north-south oriented X-Y mounts, 

^ = sin E/(l — cos 2 6 cos 2 £) 1/2 ( 4 - 47 ) 

where E is the elevation angle and 9 the azimuth (E of N). Finally, for east-west oriented X-Y mounts, 

^ = sin E/(l — sin 2 8 cos 2 E) 1 ^ 2 ( 4 - 48 ) 


4.7 OBSERVATION CORRELATIONS VIA THE TROPOSPHERE 

An alternative description of the effects of the turbulent non-hydrostatic part of the troposphere 
is obtained by introducing correlations between observations into the least-squares formalism. The 
mathematical details of this method will be described by Branson (1994). It is based on the paper of 
TVeuhaft and Lanyi (1987), as implemented by Edwards (1988). A future revision of this document 
will include a summary of the model for tropospheric correlations. For the present, Table XI gives the 
default values of the Treuhaft-Lanyi parameters that are used in MODEST. It should be emphasized 
that the present MODEST code is only capable of treating single-baseline observing sessions. 


Table XI 


Troposphere Turbulence Parameters 


Parameter 

Symbol 

Value 

Scale factor 

c 

1.7 x 10 -7 m -1 / 3 

Wind speed 

V 

8 m/s 

Wind direction 

... 

East 

Wet troposphere scale height 

h 

2 km 

Decorrelation length 

L 

3000 km 
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SECTION 5 


IONOSPHERE MODEL 


The second component of the Earth’s atmosphere, the ionosphere, is a layer of plasma at about 
350 km altitude, created primarily by the ultraviolet portion of the sunlight. In the quasi-longitudinal 
approximation (Spitzer, 1962) the refractive index of this medium is 

<>■» 

where the plasma frequency, v p , is 

/ \ 1/2 
v p = (pe 2 r 0 /ffj « 8. 

the electron gyrofrequency, v Qy is 

- tB 
U ° 2 tt me 

and © is the angle between the magnetic field B and the direction of propagation of the wave front. 
Here p is the number density of the electrons, e and m are the electron charge and mass, ro is the 
classical electron radius, and c is the speed of light. 

Tables XII and XIII give the plasma, i/ p} and gyro frequencies, u gi for the three regimes of a radio 
signal’s ray path: Earth, interplanetary, and interstellar space. Given the S-band (vs = 2.3 GHz) and 
X-band (i/x = 8-4 GHz) frequencies typically used for geodetic and astrometric VLBI, the importance 
of correcting for plasma and gyrofrequency can be parameterized by the ratios of v p and u g to and 
to ux respectively. 


.97 x 10V /2 


(5.2) 


(5.3) 


Table XII 
Plasma Effects 


Plasma 

p («/ m3 ) 

u p (kHz) 

{"p/vs) 

(vp/vx) 

Earth’s ionosphere 

10 12 

8900 

4 x 10 -3 

nr 3 

Interplanetary 

10 7 - 10 8 

28-89 

4 x 1(T 5 

10~ 5 

Interstellar 

10 5 

3 

1.2 x 10“ c 

3 x icr 7 


Table XIII 

Electron Gyrofrequency Effects 


Magnetic field 

B (gauss) 


KM) 

KM-) 

Earth 

0.2 

600 

3 x io -4 

7.5 x 10" s 

Interplanetary 

i 

o 

0.3 

1.5 x io -7 


Interstellar 

10 _G 

0.003 

1.5 x 10 -s 

3.2 x IO -10 
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Relative to vacuum as a reference, the phase delay of a monochromatic signal transiting this 
medium of refractive index n is 

3 


Tpd 


-i /(.-.>«~£/(*)Vi® ,+ i(# + -] j ' 


where 


1 - 1/2 


For 8.4 GHz, we may approximate this effect to parts in 10 c — 10 by. 




(^)cos© 


1± 


where 


,-sa / pdi- 

2* J 


crple 

2ir 


(5.4) 

(5.5) 

(5.6) 

(5.7) 


and where L is the total number of electrons per unit area along the integr ated line of sight. The 
Jar symbolise, a geometrical average. If we also neglect the term (v,cos © )Ju then the expression 
for A p d becomes simple and independent of the geometry of the traversal of the wave front throug 

the ionosphere: A pd = - 9 /^ 2 M 

This delay is negative. Thus, a phase advance actually occurs for a monochromatic signal. Since phase 
delay is obtained at a single frequency, observables derived from phase delay (c.g. phase d«Uy rates) 
experience an increment which is negative (the observable with the medium present is smaller than it 
would be without the medium). In contrast, group delays measured by a technique such as bandwidth 

synthesis (r = |^) experience an additive delay which can be derived from (5.8) by differentiating 


(5.9) 


4> = vApd with respect to frequency: 

A gd = g/l' 

Notice that the sign is now positive, though the group delay is of the same magnitude as the phase 
delay advance. For group delay measurements, the measured delay is larger with the medium present 

than without the medium. . v 

For a typical ionosphere, r * 1 to 20 xlO" 10 s at local zenith for i/ = 8.4 GHz This effect has 

a maximum at approximately 1400 hours local time and a broad minimum during local night, tor 
long baselines, the effects at each station are quite different. Thus, the differential effect may be o 

the same order as the maximum. . , 

For the interplanetary medium and at an observing frequency of 8.4 GHz, a single ray path 
experiences a delay of approximately 6x 10" 7 s in transiting the Solar System. However, the differential 
between the ray paths to the two stations on the Earth is considerably less, since the gradient between 
the two ray paths should also be inversely proportional to the dimensions of the plasma region (e.g., 
one astronomical unit as opposed to a few thousand kilometers) The ray path from a source at a 
distance of 1 megaparsec (3 x 10 7 km) experiences an integrated plasma delay of approximately 5000 
seconds for a frequency of 8.4 GHz. In this case, however, the typical dimension is also that much 
greater, and so the differential effect on two ray paths separated by one Earth radius is still not as 
great as the differential delays caused by the Earth’s ionosphere. 


5.1 


DUAL-FREQUENCY CALIBRATION 


These plasma effects can best be removed by the technique of observing the sources at two 
frequencies, v x and u 2 , where v lja » u r and where \v 2 - Vil/fa + "i) « 1- Then at the two 
frequencies v\ and u 2 we obtain 

t„i = r + <?/!/? (° 10 ) 
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and 


T v i = r + q/vl (5.11) 

Multiplying each expression by the square of the frequency involved and subtracting, we obtain 


T — G7V 2 “f" bT^i 


(5.12) 


where 


and 



(5.13) 


(5.14) 


This linear combination of the observables at two frequencies thus removes the charged particle con- 
tribution to the delay. 

For uncorrelated errors in the frequency windows, the overall error in the derived delay can be 
modeled as 


_ 2_2 


+ b 2 a\ 


2 

T* l 


(5.15) 


Modeling of other error types is more difficult and will not be treated in this report. Since the values 
of a and b are independent of q, these same coefficients apply both to group delay and to phase delay. 

If we had not neglected the effect of the electron gyrofrequency in the ionosphere, then instead 
of (5.12) above, we would have obtained 


t — ■+■ bT u \ H- 


q i/ g cos © 
Wife - ^l) 


(5.16) 


where a and b are defined as in (5.13) and (5.14), respectively. 

If we express the third term on the right-hand side in units of the contribution of the ionosphere 
at frequency u 2y we obtain 


t — ar u 2 + br^i -f- 


Aprfl/2 ^gCOS© 
Vi(l / 2 + Vl) 


(5.17) 


For X band A p< * w 1 to 20 xlO 10 s at the zenith. When using S band as the other frequency in the 
pair, this third term is % 2 x 10~ 4 A pd cos© as 2 to 40 xl(T 13 s at zenith. In the worst case of high 
ionospheric electron content, and at low elevation angles, this effect could reach 0.1 cm of total error 
in determining the total delay using the simple formula (5.12) above. Notice that the effect becomes 
much more significant at lower frequencies. 

In the software chain used at the Jet Propulsion Laboratory, the dual-frequency correction is 
performed prior to the processing step “MODEST” (Lowe, 1992). MODEST does not have the 
facility to perform this correction. However, the process is described here because it is important to 
understanding the data input to MODEST. For millimeter accuracy, or for lower observing frequencies 
even at centimeter accuracy levels, a correction for the gyrofrequency effect is necessary. 


5.2 TOTAL ELECTRON CONTENT 


In the absence of the dual-frequency observation capability described above, one can improve 
the model of the interferometer by modeling the ionosphere, using whatever measurements of the total 
electron content are available. The model we have chosen to implement is very simple. Its formalism 
is similar to that of the troposphere model, except that the ionosphere is modeled as a spherical shell 
for which the bottom is at the height above the geodetic surface of the Earth, and the top of the 
shell is at the height h 2i above that same surface (see Figure 7). For each station the ionospheric 
delay is modeled as 


t; = kgl e S(E)/v 2 


(5.18) 
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where 


k = 


O.lcrp 

27T 


(5.191 


I, b the total electron content at eenith (in electron, per meter eq«»«d xlO- 1 ’), ^ l [j°' 

Lop (phaae) delay. E la the apparent geodetic elevation angle ot the eoorce, S(E) » a alant range 

factor discussed below, and u is the observing frequency in gigaherti. 

The slant range factor (see Figure 7) is 


\/<R 2 sin 2 £ + 2flh 2 + h! - \j R 2 sin 2 E + 2Rhi + h? 

S(E) = ^ J^Thl 


(5.20) 



Figure 7. The geometry of the spherical ionospheric shell used for ionospheric corrections 

This expression is strictly correct for a spherical Earth of radius R, and a source at apparent elevation 
ancle E The model employed uses this expression, with the local rad.us of curvature of the geo.d 
surface at the receiving station R taken to be equal to the distance from the station to the center of 
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the Earth. The model also assumes this same value of R can be used at the ionospheric penetration 
points, e.g.: 

R, = R + hi (5.21) 

This is not strictly true, but is a very close approximation, particularly compared to the crude nature 
of the total electron content determinations on which the model also depends. The total ionospheric 
contribution on a given baseline is 


^ 3 1 (5.22) 

We assume that the ionospheric total electron content, J e , is the sum of two parts, one obtained by 
some external set of measurements such as Faraday rotation or GPS techniques, and the other by 
some specified additive constant: 

Ze = I' meat “h Ze add (5.23) 

These external measurements, in general, are not along directions in the ionosphere coincident with 
the ray paths to the interferometer. Thus, for each antenna, it is necessary to map a measurement 
made along one ray path to the ray paths used by the interferometer. Many different techniques to 
do this mapping have been suggested and tried; all of them are of dubious accuracy. In the light 
of these problems, and in the anticipation that dual-frequency observations will be employed for the 
most accurate interferometric work, we have implemented only a simple hour-angle mapping of the 
time history of the measurements of I e at a given latitude and longitude to the point of interest. In 
this model we allow the user to adjust the “height”, h , of the ionosphere, but require 

h\ = h — 35 km 

^2 = h + 70 km (5.24) 

Nominally, this height is taken to be 350 km. Setting this height to zero causes the program to 
ignore the ionosphere model, as is required if dual-frequency observations have already been used to 
remove the plasma effects. As in the troposphere model, these corrections can also be incorporated 
into the input data stream. Then the user is free to accept the passed correction, and use this model 
as a small alteration of the previously invoked model, or to remove the passed corrections. 

The deficiencies of these ionosphere models for single-frequency observations are compounded 
by the lens effect of the solar plasma. In effect, the Solar System is a spherical plasma lens which 
will cause the apparent positions of the radio sources to be shifted from their actual positions by an 
amount which depends on the solar weather and on the Sun-Earth-source angle. Since both the solar 
weather and the Sun-Earth-source angle change throughout the year, very accurate observations over 
the time scale of a year will be virtually impossible. 

Only one parameter is present in the ionosphere portion of the model. Again, the model is linear 
in the parameter I e add . Thus, the partial derivative with respect to this parameter is 

dr k /(station #) g { data type) 5(E) 

die add ^ (5 - 25 ) 

with /( 2) = 1 and /( 1) = —1. 
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SECTION 0 


MODELING THE PHASE DELAY RATE (FRINGE FREQUENCY) 

The interferometer is capable of producing several data types: group delay, phase delay, and the 
time rate of change of phase delay. Actually, the time rate of change of group delay is also available. 
However, it is not accurate enough to be of significance for geodetic uses. The models discussed above 
are directly applicable either to group delay or to phase delay. However, the model for the time rate 
of change of phase delay (fringe frequency) must be either constructed separately, or its equivalent 
information content obtained by forming the time difference of two phase delay values constructed 
from the delay-rate measurements as shown below. We chose the latter route since then only models of 
delay are needed. The two phase delay values, T pd {t±A), used to represent the delay-rate measurement 
information content are obtained from the expression 

r pd (t±A) = r m (t±A) + r r (t)±T r A (61) 

where r m (t) is the model used in the delay extraction processing step, r r (t) is the residual of the 
observations from that model, and f r is the residual delay rate of the data relative to that model. 
This modeling for the delay extraction step is covered in Thomas (1981), and is done in analysis 
steps prior to and completely separate from the modeling described in this report. The output of 
those previous steps is such that the details of all processing prior to the modeling described here 
are transparent to this step. Only total interferometer delays and differenced total interferometer 
phase delays (these phase delays are divided by the time interval of the difference) are reported to 
this step. One of the requirements of these previous processing steps is that the model delay used 
be accurate enough to provide a residual phase that is a linear function of time over the observation 
interval required to obtain the delay information. A linear fit to this residual phase yields the value 
of T r , the residual delay rate. Using these two values of r pd , obtained by (6.1) above, the quantity, R, 
is constructed by the following algorithm: 

R , k P d(t + A) -^r p<i (t ~ A)) ( 6 .2) 

This value and the group delay measurement, are the two data types that normally serve as 
the interferometer data input to be explained by the model described in this report. The software, 
however, also has the option to model phase delay, r p d, directly. In the limit A — ♦ 0, this expression 
for differenced phase delay approaches the instantaneous time rate of change of phase delay (fringe 
frequency) at time t . In practice, A must be large enough to avoid roundoff errors that arise from 
taking small differences of large numbers, but should also be small enough to allow R to be a reasonably 
close approximation to the instantaneous delay rate. A suitable compromise appears to be A « 2 
seconds. Fortunately, A has a wide range of allowed values, and the capability to model interferometer 
performance accurately is relatively insensitive to this choice. 
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SECTION 7 


PHYSICAL CONSTANTS USED 


In the software that has been implemented we have tried to use the constants recommended by 
the IERS Standards (1992). Those that have not been previously defined in the text, but which have 
an effect on the results that are obtained using the JPL software, are given below: 


Symbol 

Value 

Quantity 

c 

299792.458 

Velocity of light (km/s) 

r 0 

2.817938 x 10“ 15 

Classical radius of the electron (meters) 

Re 

6378.140 

Equatorial radius of the Earth (km) 

u>£ 

7.2921151467 x 10“ 6 

Rotation rate of the Earth (rad/s) 

f 

298.257 

Flattening factor of the geoid 

h 2 

0.609 

Vertical quadrupole Love number 

k 

0.0852 

Horizontal quadrupole Love number 

h 3 

0.292 

Vertical octupole Love number 

h 

0.0151 

Horizontal octupole Love number 

9 

980.665 

Surface acceleration due to gravity (cm/s 2 ) 

1PPN 

1 

Post-Newtonian relativistic gamma factor 


Another group of constants is read from the planetary ephemeris file (presently JPL DE200) when 
MODEST is executed. These constants include all planetary masses and radii; those most important 
in VLBI modeling are: 

AU 
GM S 
GM e 
Rs 

M b /M m 


1.4959787066 x 10 s 
1.3271243993544841 x 10 11 
398600.4480734463 
696000.0 
81.300587 


Astronomical unit (km) 
Mass of Sun (km 3 /s 2 ) 
Mass of Earth (km 3 /s 2 ) 
Radius of Sun (km) 
Earth/Moon mass ratio 
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SECTION 8 


POSSIBLE IMPROVEMENTS TO THE CURRENT MODEL 


This section lists areas in which the current model can be improved. 

General Relativity: 

Second-order effects have not been carefully investigated, and could possibly contribute at 
the picosecond level. 

Variations of the Earth’s gravitational potential must be taken into account m defining proper 
lengths. This correction is estimated by Thomas (1991) to amount to 0.2 cm for a 10,000 km 

baseline. 


Earth Tidal Models: 

In addition to the eight frequencies considered in the JPL92 model of Section 2.6. 1.1.2, 
short-period variations of UTPM may have components significant at the mm level that will 
emerge as data analyses are refined. 

Empirical estimates of ocean loading amplitudes for several IRIS stations (Sovers, 1994) 
indicate that the best theoretically derived amplitudes might be in error by several mm. 
Future refinements in data analyses are expected to improve the accuracy of the ocean 
loading model to the mm level. 

Resonance with the free core nutation may modify some of the amplitude corrections for 
nearly diurnal frequencies by wl mm. 

Tides cause motion of the center of mass of the solid Earth due to motion of the center of 
mass of the oceans (Brosche and Wiinsch, 1993). The amplitude of this displacement can be 
as large as 1 cm at the usual diurnal and semidiurnal tidal frequencies. Its effect on L 
observations must be assessed. 

The retarded tidal potential effect mentioned in Section 2.4.1 can be as large as several 
tenths of mm. Thus, for correct modeling at the mm level, the light travel time should be 
accounted for. 


Source Structure: 

Estimates of parameters for simplified structural models might provide improved data anal 
yses via “poor man’s mapping”. 


Earth Orientation Models: 

There are short-period deficiencies in the present IAU models for the orientation of the Earth 
in space that may be as large as 1 to 2 milliarcseconds, and longer-term deficiencies of the 
order of 1 milliarcsecond per year (3 cm at one Earth radius). VLBI measurements made 
during the past decade indicate the need for revisions of this order of the annual nutation 
terms and the precession constant [Eubanks et al. (1985), Herring et al. (1986)]. The 18.6- 
year components of the IAU nutation series are also in error, and present data spans are 
just approaching durations long enough to separate them from precession. To provide an 
improved nutation model, we have implemented MODEST options to use the nutation series 
of Zhu et al., ZMOA 1990-2, or Kinoshita and Souchay (discussed in Section 2.6. 2.1). Any 
of these constitute a provisionally improved model, especially for the annual and semiannual 
nutations, until the IAU series is officially revised. 
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Antenna Deformation: 


Gravity loading may cause variations in the position of the reference point of a large antenna 
that are as large as 1 cm in the local vertical direction. Liewer (1986) presents evidence 
that this causes systematic errors and that their dependence on antenna orientation and 
temperature may be modeled. 


Antenna Alignment: 

Hour angle misalignment on the order of 1 arc minute can cause 1 mm delay effects for DSN 
HA-Dec antennas with 7-m axis offsets. 


Troposphere: 

New techniques for characterizing the atmosphere are expected to allow more realistic model- 
ing of the tropospheric delay than the simple spherical-shell model underlying all the results 
of Section 4. When comprehensive atmospheric data from a region surrounding each observ- 
ing site are available, present computer speeds permit estimating the tropospheric delay by 
means of a complete ray-tracing solution. Meanwhile, improvements in tropospheric map- 
ping should be sought by modeling details such as the global, seasonal, and diurnal variation 
of the inversion layer height and troposphere azimuthal asymmetry. 


Ionosphere: 

Corrections for the gyrofrequency effect may be of millimeter order. 

Thermal Effects: 

Thermal expansion of the portion of an antenna above the reference point will induce delay 
signatures that are « 9 ps (3 mm) peak-to-peak for a 34-m dish. 


Phase Delay Rate: 

Rather than modeling the delay rates as finite differences of model delays, direct analytic 
expressions for derivatives of delays could be implemented. This would eliminate questions 
concerning the choice of the time difference A discussed in Section 6. Care must be exercised, 
however, to ensure consistency between definitions of modeled and observed delay rates. 
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appendix a 


NUTATION MODELS 

Tb. five nutntion ..tie. nvnilnbl. in MODEST nr. tnbnlnt.fi h.r.t T.bl. A I give. th. .t.nfiyfi 1980 
IAU series- Tables A.II, A.III, and A.IV contain the results of Zhu et al. (1990), Tables A. V, A. , 
nnd 1 VH contnin th, ZMOA 1990-2 T.bl.t A. VIII nnfi A.DC cnntn.n the K.no.h.tn-Snnohn, 

series; for completeness, the old (Woolard) nutation series is given in Table A.X. 

Table A.I 

1980 IAU Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 
kji kf 2 kj 3 kji fcys 

Aoj Al, 

( 0 ". 0001 ) 
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46 
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3232.9 
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48 
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Table A. I cont. 


1980 IAU Theory of Nutation 


Index 

j 

Period 

(days) 

k 3l 

Argument coefficient 
kj 2 fc/3 kj 4 


Aoj Aij 
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Boy B\j 
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31 
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2 

0 

2 

-31 

0.0 

13 

0.0 

46 

13.6 

0 

0 

2 

0 

0 

26 

0.0 

-1 

0.0 

47 

27.0 

-1 

0 

2 

0 

1 

21 


-10 

0.0 

48 

32.0 

-1 

0 

0 

2 


16 

0.0 

-8 

0.0 

49 

31.7 

1 

0 

0 

-2 


-13 

0.0 

7 

0.0 

50 

9.5 

-1 

0 

2 

2 


-10 

0.0 

5 

mm 

51 

34.8 

1 

1 

0 

-2 


-7 

0.0 

0 

0.0 

52 

13.2 

0 

1 

2 

0 

2 

7 

0.0 

-3 

0.0 

53 

14.2 

0 

-1 

2 

0 

2 

-7 

0.0 

3 

0.0 

54 

5.6 

1 

0 

2 

2 

2 

-8 

0.0 

3 

0.0 

55 

9.6 

1 

0 

0 

2 

0 

6 

0.0 

0 

0.0 

56 

12.8 

2 

0 

2 

-2 

2 

6 

0.0 

-3 

0.0 

57 

14.8 

0 

0 

0 

2 

1 

-6 

0.0 

3 

0.0 

58 

7.1 

0 

0 

2 

2 

1 

-7 

0.0 

3 

0.0 

59 

23.9 

1 

0 

2 

-2 

1 

6 

0.0 

-3 

0.0 

60 

14.7 

0 

0 

0 

-2 

1 

-5 

0.0 

3 

0.0 

61 

29.8 

1 

-1 

0 

0 

0 

5 

0.0 

0 

0.0 

62 

6.9 

2 

0 

2 

0 

1 

-5 

0.0 

3 

0.0 

63 

15.4 

0 

1 

0 

-2 

0 

-4 

0.0 

0 

0.0 

64 

26.9 

1 

0 

-2 

0 

0 

4 

0.0 

0 

0.0 

65 

29.5 

0 

0 

0 

1 

0 

-4 

0.0 

0 

0.0 

66 

25.6 

1 

1 

0 

0 

0 

-3 

0.0 

0 

0.0 

67 

9.1 

1 

0 

2 

0 

0 

3 

0.0 

0 

0.0 

68 

9.4 

1 

-1 

2 

0 

2 

-3 

0.0 

1 

0.0 

69 

9.8 

-1 

-1 

2 

2 

2 

-3 

0.0 

1 

0.0 

70 

13.7 

-2 

0 

0 

0 

1 

-2 

0.0 

1 

0.0 


80 










Table A . I cont . 


1980 IAU Theory of Nutation 


Index 

Period 

Argument coefficient 


Aoj 

Mi 

Boy 

Bxi 

j 

( days ) 

fcyi 

kj2 

kjz 


fcj'5 

( 0 ". 0001 ) 

(V .UUU1J 

71 

72 

73 

74 

75 

5.5 

7.2 

8.9 

32.6 

13.8 

3 

0 

1 

-1 

2 

0 

-1 

1 

0 

0 

2 

2 

2 

2 

0 

0 

2 

0 

-2 

0 

2 

2 

2 

1 

1 

-3 

-3 

2 

-2 

2 

0.0 

0.0 

0.0 

0.0 

0.0 

1 

1 

-1 

1 

-1 

0.0 

0.0 

0.0 

0.0 

0.0 

76 

27.8 

1 

0 

0 

0 

2 

-2 

0.0 

1 

0.0 

77 

9.2 

3 

0 

0 

0 

0 

2 

0.0 

0 

0.0 

78 

9.3 

0 

0 

2 

1 

2 

2 

0.0 

-1 

0.0 

79 

27.3 

-1 

0 

0 

0 

2 

1 

0.0 

-1 

0.0 

80 

10.1 

1 

0 

0 

-4 

0 

-1 

0.0 

0 

0.0 

81 

14.6 

-2 

0 

2 

2 

2 

1 

0.0 

•1 

0.0 

82 

5.8 

-1 

0 

2 

4 

2 

-2 

0.0 

1 

0.0 

83 

15.9 

2 

0 

0 

-4 

0 

-1 

0.0 

0 

0.0 

84 

22.5 

1 

1 

2 

-2 

2 

1 

0.0 

-1 

0.0 

85 

5.6 

1 

0 

2 

2 

1 

-1 

0.0 

1 

0.0 

86 

7.3 

-2 

0 

2 

4 

2 

-1 

0.0 

1 

0.0 

87 

9.1 

-1 

0 

4 

0 

2 

1 

0.0 

0 

0.0 

88 

29.3 

1 

-1 

0 

-2 

0 

1 

0.0 

0 

0.0 

89 

12.8 1 

2 

0 

2 

-2 

1 

1 

0.0 

.1 

0.0 

90 

4.7 

2 

0 

2 

2 

2 

-1 

0.0 

0 

0.0 

91 

9.6 

1 

0 

0 

2 

1 

-1 

0.0 

0 

0.0 

92 

12.7 

0 

0 

4 

-2 

2 

1 

0.0 

0 

0.0 

93 

8.7 

3 

0 

2 

-2 

2 

1 

0.0 

0 

0.0 

94 

23.8 

1 

0 

2 

-2 

0 

-1 

0.0 

0 

0.0 

95 

13.1 

0 

1 

2 

0 

1 

1 

0.0 

0 

0,0 

96 

35.0 

-1 

-1 

0 

2 

1 

1 

0.0 

0 

0.0 

97 

13.6 

0 

0 

-2 

0 

1 

-1 

0.0 

0 

0.0 

98 

25.4 

0 

0 

2 

-1 

2 

-1 

0.0 

0 

0.0 

99 

14.2 

0 

1 

0 

2 

0 

-1 

0.0 

0 

0.0 

100 

9.5 

1 

0 

-2 

-2 

0 

-1 

0.0 

0 

0.0 

101 

14.2 

0 

-1 

2 

0 

1 

-1 

0.0 

0 

0.0 

102 

34.7 

1 

1 

0 

-2 

1 

-1 

0.0 

0 

0.0 

103 

32.8 

1 

0 

-2 

2 

0 

-1 

0.0 

0 

0.0 

104 

7.1 

2 

0 

0 

2 

0 

1 

0.0 

0 

0.0 

105 

4.8 

0 

0 

2 

4 

2 

-1 

0.0 

0 

0.0 

106 

27.3 

0 

i 

1 

0 

1 

0 

1 

j 

0.0 

0 

0.0 


81 



Table A. II 


Zhu et at. Theory of Nutation: 1980 IAU Terms 


Index 

j 

Period 

(days) 

kji 

Argument coefficient 
kj 2 k :4 k 

ys 

Aoy Ai j 

(0".00001) 

Bo: B\j 

(0",00001) 

1 

6798.38 

0 

0 

0 

0 

l 

-1720618 

-1743 

920530 

90 

2 

3399.19 

0 

0 

0 

0 

2 

20743 

2 

-8975 

5 

3 

1305.48 

-2 

0 

2 

0 

l 

460 

1 

-243 

0 

4 

1095.18 

2 

0 

-2 

0 

0 

110 

0 

1 

0 

5 

1615.75 

-2 

0 

2 

0 

2 

-31 

0 

14 

0 

6 

3232.86 

1 

-1 

0 

-1 

0 

-33 

0 

0 

0 

7 

6786.32 

0 

-2 

2 

-2 

1 

-15 

0 

8 

0 

8 

943.23 

2 

0 

-2 

0 

1 

7 

0 

-4 

0 

9 

182.62 

0 

0 

2 

-2 

2 

-131720 

-16 

57320 

-31 

10 

365.26 

0 

1 

0 

0 

0 

14735 

-35 

719 

-2 

11 

121.75 

0 

1 

2 

-2 

2 

-5176 

12 

2247 

-7 

12 

365.22 

0 

-1 

2 

-2 

2 

2161 

-5 

-961 

3 

13 

177.84 

0 

0 

2 

-2 

1 

1293 

1 

-699 

0 

14 

205.89 

2 

0 

0 

-2 

0 

479 

0 

5 

0 

15 

173.31 

0 

0 

2 

-2 

0 

-218 

0 

-1 

0 

16 

182.63 

0 

2 

0 

0 

0 

168 

-1 

2 

0 

17 

386.00 

0 

1 

0 

0 

1 

-140 

0 

86 

0 

18 

91.31 

0 

2 

2 

-2 

2 

-158 

1 

69 

0 

19 

346.64 

0 

-1 

0 

0 

1 

-127 

0 

64 

0 

20 

199.84 

-2 

0 

0 

2 

1 

-58 

0 

30 

0 

21 

346.60 

0 

-1 

2 

-2 

1 

-48 

0 

27 

0 

22 

212.32 

2 

0 

0 

-2 

1 

41 

0 

-22 

0 

23 

119.61 

0 

1 

2 

-2 

1 

36 

0 

-20 

0 

24 

411.78 

1 

0 

0 

-1 

0 

-43 

0 

-6 

0 

25 

131.67 

2 

1 

0 

-2 

0 

11 

0 

0 

0 

26 

169.00 

0 

0 

-2 

2 

1 

9 

0 

-4 

0 

27 

329.79 

0 

1 

-2 

2 

0 

-9 

0 

0 

0 

28 

409.23 

0 

1 

0 

0 

2 

7 

0 

-3 

0 

29 

388.27 

-1 

0 

0 

1 

1 

9 

0 

-4 

0 

30 

117.54 

0 

1 

2 

-2 

0 

-6 

0 

0 

0 

31 

13.66 

0 

0 

2 

0 

2 

-22824 

-2 

9806 

-5 

32 

27.55 

1 

0 

0 

0 

0 

7122 

1 

-70 

0 

33 

13.63 

0 

0 

2 

0 

1 

-3885 

-4 

2011 

0 

34 

9.13 

1 

0 

2 

0 

2 

-3023 

0 

1293 

-1 

35 

31.81 

1 

0 

0 

-2 

0 

-1572 

0 

-13 

0 

36 

27.09 

-1 

0 

2 

0 

2 

1238 

0 

-535 

0 

37 

14.77 

0 

0 

0 

2 

0 

635 

0 

-13 

0 

38 

27.67 

1 

0 

0 

0 

1 

633 

1 

-332 

0 

39 

27.44 

-1 

0 

0 

0 

1 

-580 

-1 

315 

0 

40 

9.56 

-1 

0 

2 

2 

2 

-598 

0 

256 

0 


82 






Table A. II cont. 


Zhu ct al. Theory of Nutation: 1980 1AU Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kji kj2 kj2 kj4 

kjB 

■^Oj ^lj 

( 0 ". 00001 ) 

Bo) B\j 

( 0 " . 00001 ) 

41 

9.12 

1 

0 

2 

0 

1 

-517 

0 

265 

0 

42 

7.10 

0 

0 

2 

2 

2 

-386 

0 

165 

0 

43 

13.78 

2 

0 

0 

0 

0 

293 

0 

-6 

0 

44 

23.94 

1 

0 

2 

-2 

2 

286 

0 

-124 

0 

45 

6.86 

2 

0 

2 

0 

2 

-311 

0 

132 

0 

46 

13.61 

0 

0 

2 

0 

0 

259 

0 

-5 

0 

47 

26.98 

-1 

0 

2 

0 

1 

205 

0 

-107 

0 

48 

31.96 

-1 

0 

0 

2 

1 

152 

0 

o 

GO 

1 

0 

49 

31.66 

1 

0 

0 

-2 

1 

-129 

0 

70 

0 

50 

9.54 

-1 

0 

2 

2 

1 

-103 

0 

53 

0 

51 

34.85 

1 

1 

0 

-2 

0 

-74 

0 

-1 

0 

52 

13.17 

0 

1 

2 

0 

2 

76 

0 

-33 

0 

53 

14.19 

0 

-1 

2 

0 

2 

-71 

0 

31 

0 

54 

5.64 

1 

0 

2 

2 

2 

-77 

0 

32 

0 1 

55 

9.61 

1 

0 

0 

2 

0 

66 

0 

-3 

0 

56 

12.81 

2 

0 

2 

-2 

2 

65 

0 

-28 

0 

57 

14.80 

0 

0 

0 

2 

1 

-64 

0 

33 

0 

58 

7.09 

0 

0 

2 

2 

1 

-66 

0 

34 

0 

59 

23.86 

1 

0 

2 

-2 

1 

58 

0 

-30 

0 

60 

14.73 

0 

0 

0 

-2 

1 

-50 

0 

28 

0 

61 

29.80 

1 

-1 

0 

0 

0 

47 

0 

-1 

0 

62 

6.85 

2 

0 

2 

0 

1 

-53 

0 

27 

0 

63 

15.39 

0 

1 

0 

-2 

0 

-44 

0 

-1 

0 

64 

26.88 

1 

0 

-2 

0 

0 

41 

0 

1 

0 

65 

29.53 

0 

0 

0 

1 

0 

-40 

0 

1 

0 

66 

25.62 

1 

1 

0 

0 

0 

-34 

0 

1 

0 

67 

9.11 

1 

0 

2 

0 

0 

34 

0 

-1 

0 

68 

9.37 

1 

-1 

2 

0 

2 

-29 

0 

12 

0 

69 

9.81 

-1 

-1 

2 

2 

2 

-29 

0 

12 

0 

70 

13.75 

-2 

0 

0 

0 

1 

* 

to 

00 

0 

13 

0 

71 

5.49 

3 

0 

2 

0 

2 

-29 

0 

12 

0 

72 

7.24 

0 

-1 

2 

2 

2 

-26 

0 

11 

0 

73 

8.91 

1 

1 

2 

0 

2 

25 

0 

-10 

0 

74 

32.61 

-1 

0 

2 

-2 

1 

-20 

0 

11 

0 

75 

13.81 

2 

0 

0 

0 

1 

22 

0 

-11 

0 

76 

27.78 

1 

0 

0 

0 

2 

-20 

0 

8 

0 

77 

9.18 

3 

0 

0 

0 

0 

16 

0 

-1 

0 

78 

9.34 

0 

0 

2 

1 

2 

16 

0 

-7 

0 

79 

27.33 

-1 

0 

0 

0 

2 

14 

0 

-6 

0 

80 

10.08 

1 

0 

0 

-4 

0 

-14 

0 

-1 

0 


83 



Table A. II cont. 


Zhu ct al. Theory of Nutation: 1980 IAU Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kji ^j‘3 


Aoj Ai, 

(0" .00001) 

•Boy Bij 

(0".0000l) 

81 

14.63 

-2 

0 

2 

2 

2 

13 

0 

-6 

0 

82 

5.80 

-1 

0 

2 

4 

2 

-15 

0 

6 

0 

83 

15.91 

2 

0 

0 

-4 

0 

-13 

0 

0 

0 

84 

22.47 

1 

1 

2 

-2 

2 

13 

0 

-5 

0 

85 

5.64 

1 

0 

2 

2 

1 

-13 

0 

7 

0 

86 

7.35 

-2 

0 

2 

4 

2 

-12 

0 

5 

0 

87 

9.06 

-1 

0 

4 

0 

2 

11 

0 

-5 

0 . 

88 

29.26 

1 

-1 

0 

-2 

0 

9 

0 

0 

0 

89 

12.79 

2 

0 

2 

-2 

1 

10 

0 

-5 

0 

90 

4.68 

2 

0 

2 

2 

2 

-11 

0 

5 

0 

91 

9.63 

1 

0 

0 

2 

1 

-10 

0 

5 

0 

92 

12.66 

0 

0 

4 

-2 

2 

9 

0 

-4 

0 

93 

8.75 

3 

0 

2 

-2 

2 

9 

0 

-4 

0 

94 

23.77 

1 

0 

2 

-2 

0 

-7 

0 

0 

0 

95 

13.41 

0 

1 

2 

0 

1 

8 

0 

-4 

0 

96 

35.03 

-1 

-1 

0 

2 

1 

7 

0 

-4 

0 

97 

13.58 

0 

0 

-2 

0 

1 

-6 

0 

3 

0 

98 

25.42 

0 

0 

2 

-1 

2 

-7 

0 

3 

0 

99 

14.19 

0 

1 

0 

2 

0 

-6 

0 

0 

0 

100 

9.53 

1 

0 

-2 

-2 

0 

-6 

0 

0 

0 

101 

14.16 

0 

-1 

2 

0 

1 

-7 

0 

3 

0 

102 

34.67 

1 

1 

0 

-2 

1 

-6 

0 

3 

0 

103 

32.76 

1 

0 

-2 

2 

0 

-6 

0 

0 

0 

104 

7.13 

2 

0 

0 

2 

0 

6 

0 

0 

0 

105 

4.79 

0 

0 

2 

4 

2 

-7 

0 

3 

0 

106 

27.32 

0 

1 

0 

1 

0 

5 

0 

0 

0 


Table A. Ill 


Zhu ct at. Theory of Nutation: Out-of-Phase Terms 


Index 

Period 


Argument coefficient 



B 2 / 

j 

(days) 




fc,4 

*ys 

(0".00001) 

1 

6798.38 

0 

0 

0 

0 

1 

221 

112 

2 

182.62 

0 

0 

2 

-2 

2 

-153 

-61 

3 

365.26 

0 

1 

0 

0 

0 

-55 

22 

4 

13.66 

0 

0 

2 

0 

2 

-5 

-2 


84 







Table A. IV 


Zhu et al. Theory of Nutation: Additional Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kjl kj 2 kj*4 kys 

Aoj Bo j 

(0".0000l) 

1 

5.49 

3 

0 

2 

0 

1 

-5 

2 

2 

5.73 

1 

-1 

2 

2 

2 

-6 

2 

3 

6.96 

0 

1 

2 

2 

2 

5 

-2 

4 

6.99 

2 

-1 

2 

0 

2 

-5 

2 

5 

7.38 

0 

0 

0 

4 

0 

5 

0 

6 

9.31 

•1 

1 

2 

2 

2 

6 

-2 

7 

9.80 

-1 

.1 

2 

2 

1 

-5 

2 

8 

9.87 

1 

-1 

0 

2 

0 

5 

0 

9 

14.83 

0 

0 

0 

2 

2 

-5 

2 

10 

29.93 

1 

-1 

0 

0 

1 

5 

-3 

11 

73.05 

0 

3 

2 

-2 

2 

-5 

2 

12 

177.84 

0 

0 

2 

-2 

1 

-9 

7 

13 

187.66 

0 

0 

2 

-2 

3 

13 

-2 

14 

3230.13 

-1 

-1 

2 

-1 

2 

13 

-5 

15 

3231.50 

-1 

0 

1 

0 

1 

15 

3 

16 

6164.10 

.1 

1 

0 



7 

-4 

17 

4.00 

3 

0 

2 



-1 

1 

18 

4.08 

1 

0 

2 




1 

19 

4.58 

4 

0 

2 

0 


-3 

1 

20 

4.68 

2 

0 

2 

2 

1 

-2 

1 

21 

4.79 

0 

0 

2 

4 

1 

-1 

1 

22 

5.56 

1 

1 



2 

1 

-1 

23 

5.80 

-1 

0 



1 

-3 

1 

24 

5.90 

-1 

-1 



2 

-2 

1 

25 

6.73 

2 

1 

2 

0 

2 

4 

-2 

26 

6.82 

0 

0 

4 

0 

2 

2 

-1 

27 

6.85 

2 

0 

2 

0 

0 

3 

0 

28 

6.98 

1 

0 

2 

1 

2 

3 

-1 

29 

7.08 

0 

0 

2 

2 

0 

4 

0 

30 

7.13 

2 

0 

0 

2 

1 

-1 

1 

31 

7.23 

0 

-1 

2 

2 

1 

-4 

2 

32 

7.34 

-2 

0 

2 

4 

1 

-2 

1 

33 

7.38 

0 

-2 

2 

2 

2 

-1 

1 

34 

7.39 

0 

0 

0 

4 

1 

-2 

1 

35 

8.68 

1 

0 

4 

-2 

2 

2 

-1 

36 

8.73 

3 

0 

2 

-2 

1 

2 

-1 

37 

8.90 

1 

1 

2 

0 

1 

4 

-2 

38 

9.05 

-1 

0 

4 

0 

1 

2 

-1 

39 

9.11 

0 

1 

2 

1 

2 

-2 

1 

40 

9.17 

-3 

0 

0 

0 

1 

-1 

1 
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Table A. IV cont. 


Zhu et al. Theory of Nutation: Additional Terms 


Index 

j 

Period 

(days) 

Argument coefficient 

kjl 


-4o; Boy 

(0".00001) 

41 

9.33 

0 

0 

2 

1 

1 

3 

-1 

42 

9.35 

1 

-1 

2 

0 

1 

-4 

2 

43 

9.60 

-1 

0 

0 

-2 

1 

-4 

3 

44 

10.07 

1 

0 

0 

-4 

1 

-1 

1 

45 

10.10 

-1 

0 

0 

4 

1 

-2 

1 

46 

10.37 

-1 

-1 

0 

4 

0 

1 

0 

47 

12.38 

2 

1 

2 

-2 

2 

3 

-1 

48 

12.64 

0 

0 

4 

-2 

1 

2 

-1 

49 

13.22 

1 

0 

2 

-1 

2 

-3 

1 

50 

13.28 

2 

1 

0 

0 

0 

-3 

0 

51 

13.63 

0 

0 

2 

0 

1 

-1 

0 

52 

13.69 

0 

0 

2 

0 

3 

2 

0 

53 

14.22 

0 

1 

0 

2 

1 

2 

-1 

54 

14.25 

1 

0 

0 

1 

0 

-3 

0 

55 

14.32 

2 

-1 

0 

0 

0 

4 

0 

56 

14.60 

-2 

0 

2 

2 

1 

2 

-1 

57 

14.70 

0 

0 

0 

-2 

2 

1 

-1 

58 

15.35 

0 

1 

0 

-2 

1 

-3 

2 

59 

15.42 

0 

-1 

0 

2 

1 

-2 

1 

60 

15.87 

2 

0 

0 

-4 

1 

-1 

1 

61 

15.94 

-2 

0 

0 

4 

1 

1 

-1 

62 

16.06 

0 

-2 

0 

2 

0 

2 

0 

63 

16.10 

0 

0 

2 

-4 

1 

-1 

1 

64 

22.40 

1 

1 

2 

-2 

1 

3 

-1 

65 

25.22 

-1 

1 

2 

0 

2 

4 

-2 

66 

25.53 

-1 

-1 

0 

0 

1 

2 

-1 

67 

25.72 

1 

1 

0 

0 

1 

-3 

2 

68 

26.77 

1 

0 

-2 

0 

1 

3 

-1 

69 

27.32 

0 

0 

1 

0 

1 

-2 

0 

70 

29.26 

-1 

-1 

2 

0 

2 

-2 

1 

71 

29.39 

-1 

1 

0 

2 

1 


1 

72 

29.40 

0 

0 

0 

-1 

1 



73 

29.66 

0 

0 

0 

1 

1 



74 

29.67 

-1 

1 

0 

• 0 

1 



75 

31.52 

1 

0 

0 

-2 

2 


-1 

76 

32.11 

-1 

0 

0 

2 

2 


2 

77 

32.45 

-1 

0 

2 

-2 

2 

• 3 

-1 

78 

35.80 

-1 

1 

2 

-2 

1 

-1 

1 

79 

38.52 

-1 

-2 

0 

2 

0 

3 

0 

80 

38.74 

1 

0 



1 

-4 

2 
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Table A. IV cont. 


Zhu et al. Theory of Nutation: Additional Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kyi kjs kj4 

fc/s 

•Ao; Boy 

(0".00001) 

121 

13.83 

2 

0 

0 

0 

2 

-1 

0 

122 

14.13 

-1 

0 

2 

1 

2 

1 

0 

123 

14.16 

0 

1 

0 

2 

-1 

-1 

0 

124 

14.76 

0 

-2 

2 

0 

2 

-1 

0 

125 

14.93 

2 

0 

-2 

2 

-1 

1 

0 

126 

15.24 

-2 

-1 

2 

2 

2 

1 

0 

127 

15.31 

-1 

0 

0 

3 

0 

-1 

0 

128 

16.63 

-2 

-1 

0 

4 

0 

1 

0 

129 

23.43 

-1 

0 

4 

-2 

2 

*1 

0 

130 

23.94 

1 

2 

0 

0 

0 

-1 

0 

131 

25.13 

-1 

1 

2 

0 

-1 

1 

0 

132 

25.32 

0 

0 

2 

-1 

1 

-1 

0 

133 

25.52 

1 

-1 

2 

-2 

1 

-1 

0 

134 

25.62 

1 

-1 

2 

-2 

2 

-1 

0 

135 

25.83 

2 

0 

0 

-1 

0 

1 

0 

136 

27.09 

-1 

2 

0 

2 

0 

-1 

0 

137 

27.32 

0 

-1 

2 

-1 

2 

1 

0 

138 

28.15 

3 

0 

-2 

0 

-1 

1 

0 

139 

29.14 

-1 

-1 

2 

0 

1 

-1 

0 

140 

29.14 

-1 

1 

0 

2 

-1 

-1 

0 

141 

31.06 

-3 

0 

2 

2 

1 

1 

0 

142 

32.45 

1 

-2 

0 

0 

0 

1 

0 

143 

34.48 

-2 

0 

0 

3 

0 

-1 

0 

144 

37.62 

-3 

0 

0 

4 

0 

1 

0 

145 

38.52 

-1 

0 

-2 

4 

-2 

-1 

0 

146 

38.96 

-1 

0 

-2 

4 

0 

-1 

0 

147 

43.06 

-1 

-1 

-2 

4 

-2 

1 

0 

148 

43.34 

1 

1 

2 

-4 

1 

-1 

0 

149 

90.10 

0 

2 

2 

-2 

1 

1 

0 

150 

96.78 

2 

0 

2 

-4 

2 

-1 

0 

151 

134.27 

2 

1 

0 

•2 

1 

1 

0 

152 

156.52 

-2 

0 

4 

-2 

2 

-1 

0 

153 

164.08 

-2 

2 

2 

0 

2 

-1 

0 

154 

187.67 

0 

2 

0 

0 

1 

-1 

0 

155 

193.56 

1 

-1 

2 

-3 

2 

1 

0 

156 

235.96 

-4 

0 

2 

2 

2 

-1 

0 
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Table A.V 


ZMOA 1990-2 Theory of Nutation: 1980 LAU Terms 


■ 

Period 

(days) 

Argument coefficient 
kj\ kj 2 fcj‘3 


Aoj Aij 

(0" .00001) 

m 

1 

6798.38 

0 

0 

0 

0 

1 

-1720670 

-1743 

920525 

90 

2 

3399.19 

0 

0 

0 

0 

2 

20751 

2 

-8977 

5 

3 

1305.48 

-2 

0 

2 

0 

1 

460 

0 

-243 

0 

4 

1095.18 

2 

0 

-2 

0 

0 

110 

0 

1 

0 

5 

1615.75 

-2 

0 

2 

0 

2 

-31 

0 

14 

0 

6 

3232.86 

1 

-1 

0 

-1 

0 

-33 

0 

0 

0 

7 

6786.32 

0 

-2 

2 

-2 

1 

-15 

0 

8 

0 

8 

943.23 

2 

0 

-2 

0 

1 

7 

0 

-4 

0 

9 

182.62 

0 

0 

2 

.2 

2 

-131714 

-16 

57305 

-31 

10 

365.26 

0 

1 

0 

0 

0 

14764 

-36 

732 

-2 

11 

121.75 

0 

1 

2 

-2 

2 

-5168 

12 

2244 

-7 

12 

365.22 

0 

-1 

2 

-2 

2 

2161 

-5 

-960 

3 

13 

177.84 

0 

0 

2 

-2 

1 

1291 

1 

-698 

0 

14 

205.89 

2 

0 

0 

-2 

0 

479 

0 

5 

0 

15 

173.31 

0 

0 

2 

-2 

0 

-219 

0 

-1 

0 

16 

182.63 

0 

2 

0 

0 

0 

168 

-1 

2 

0 

17 

386.00 

0 

1 

0 

0 

1 

-140 

0 

86 

0 

18 

91.31 

0 

2 

2 

-2 

2 

-158 

i 

69 

0 

19 

346.64 

0 

-1 

0 

0 

1 

-128 

0 

64 

0 

20 

199.84 

-2 

0 

0 

2 

1 

-58 

0 

31 

0 

21 

346.60 

0 

-1 

2 

-2 

1 

-48 

0 

28 

0 

22 

212.32 

2 

0 

0 

-2 

1 

41 

0 

-22 

0 

23 

119.61 

0 

1 

2 

-2 

1 

36 

0 

-20 

0 

24 

411.78 

1 

0 

0 

-1 

0 

-45 

0 

-7 

0 

25 

131.67 

2 

1 

0 

-2 

0 

11 

0 

0 

0 

26 

169.00 

0 

0 

-2 

2 

1 

9 

0 

-4 

0 

27 

329.79 

0 

1 

-2 

2 

0 

-9 

0 

0 

0 

28 

409.23 

0 

1 

0 

0 

2 

7 

0 

-3 

0 

29 

388.27 

-1 

0 

0 

1 

1 

10 

0 

-4 

0 

30 

117.54 

0 

1 

2 

-2 

0 

-6 

0 

0 

0 

31 

13.66 

0 

0 

2 

0 

2 

-22756 

-2 

9779 

-5 

32 

27.55 

1 

0 

0 

0 

0 

7111 

1 

-67 

0 

33 

13.63 

0 

0 

2 

0 

1 

-3871 

-4 

2006 

0 

34 

9.13 

1 

0 

2 

0 

2 

-3012 

0 

1289 

-1 

35 

31.81 

1 

0 

0 

-2 

0 

-1570 

0 

-13 

0 

36 

27.09 

-1 

0 

2 

0 

2 

1235 

0 

-533 

0 

37 

14.77 

0 

0 

0 

2 

0 

636 

0 

-12 

0 

38 

27.67 

1 

0 

0 

0 

1 

632 

1 

-332 

0 

39 

27.44 

-1 

0 

0 

0 

1 

-581 

-1 

315 

0 

40 

1 

9.56 

-1 

J 

0 

2 

2 

2 

-596 

0 

255 

0 
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Table A.V cont. 


ZMOA 1990-2 Theory of Nutation: 1980 IAU Terms 


Index 

j 

Period 

(days) 

kji 

Argument coefficient 

^>2 ^j‘3 

fcys 

■^Oj ^l; 

(0".00001) 

Boj 

(0".00001) 

41 

9.12 


0 

2 

0 

1 

-515 

0 

264 

0 

42 

7.10 



2 

2 

2 

-384 

0 

164 

0 

43 

13.78 



0 

0 

0 

293 

0 

-6 

0 

44 

23.94 

1 

0 

2 

-2 

2 

286 

0 

-124 

0 

45 

6.86 

2 

0 

2 

0 

2 

-309 

0 

132 

0 

46 

13.61 

0 

0 

2 

0 

0 

259 

0 

-5 

0 

47 

26.98 

-1 

0 

2 

0 

1 

205 

0 

-107 

0 

48 

31.96 

-1 

0 

0 

2 

1 

152 

0 

-80 

0 

49 

31.66 

1 

0 

0 

-2 

1 

-129 

0 

70 

0 

50 

9.54 

-1 

0 

2 

2 

1 

-102 

0 

52 

0 

51 

34.85 

1 

1 

0 

-2 

0 

-74 

0 

-1 

0 

52 

13.17 

0 

1 

2 

0 

2 

76 

0 

-33 

0 

53 

14.19 

0 

-1 

2 

0 

2 

-76 

0 

33 

0 

54 

5.64 

1 

0 

2 

2 

2 

-77 

0 

32 

0 

55 

9.61 

1 

0 

0 

2 

0 

66 

0 

-2 

0 

56 

12.81 

2 

0 

2 

-2 

2 

65 

0 

-28 

0 

57 

14.80 

0 

0 

0 

2 

1 

-64 

0 

33 

0 

58 

7.09 

0 

0 

2 

2 

1 

-66 

0 

34 

0 

59 

23.86 

1 

0 

2 

-2 

1 

58 

0 

-30 

0 

60 

14.73 

0 

0 

0 

-2 

1 

-50 

0 

28 

0 

61 

29.80 

1 

-1 

0 

0 

0 

47 

0 

-1 

0 


6.85 

2 

0 

2 

0 

1 

-53 

0 

26 

0 

■a 

15.39 

0 

1 

0 

-2 

0 

-44 

0 

-1 

0 

64 

26.88 

1 

0 

-2 

0 

0 

41 

0 

1 

0 

65 

29.53 

0 

0 

0 

1 

0 


0 

1 

0 

66 

25.62 

1 

1 

0 

0 ’ 

0 

-34 

0 

1 

0 

67 

9.11 

1 

0 

2 

0 

0 

34 

0 

-1 

0 

68 

9.37 

1 

-1 

2 

0 

2 

-29 

0 

12 

0 

69 

9.81 

-1 

-1 

2 

2 

2 

-29 

0 

12 

0 

70 

13.75 

-2 

0 

0 

0 

1 

-23 

0 

13 

0 

71 

5.49 

3 

0 

2 

0 

2 

-29 

0 

12 

0 

72 

7.24 

0 

-1 

2 

2 

2 

-26 

0 

11 

0 

73 

8.91 

1 

1 

2 

0 

2 

24 

0 

-10 

0 

74 

32.61 

-1 

0 

2 

-2 

1 

-20 

0 

11 

0 

75 

13.81 

2 

0 

0 

0 

1 

21 

0 

-11 

0 

76 

27.78 

1 

0 

0 

0 

2 

-20 

0 

8 

0 

77 

9.18 

3 

0 

0 

0 

0 

16 

0 

-1 

0 

78 

9.34 

0 

0 

2 

1 

2 

16 

0 

-7 

0 

79 

27.33 

-1 

0 

0 

0 

2 

14 

0 

-6 

0 

80 

10.08 

1 

0 

0 

-4 

0 

-14 

0 

-1 

0 
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Table A.V cont. 


ZMOA 1990-2 Theory of Nutation: 1980 IAU Terms 
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Table A. VI 


ZMOA 1990-2 Theory of Nutation: Additional Terms 


Index 

j 


Argument coefficient 

kji kj2 kjd kj 4 

kjS 

Aoy Boy 

(0".00001) 

I 

5.49 

3 

0 

2 

0 

1 

-5 

2 

2 

5.73 

1 

-1 

2 

2 

2 

-6 

2 

3 

6.96 

0 

1 

2 

2 

2 

5 

-2 

4 

6.99 

2 

-1 

2 

0 

2 

-5 

2 

5 

7.38 

0 

0 

0 

4 

0 

5 

0 

6 

9.31 

-1 

1 

2 

2 

2 

6 

-2 

7 

9.80 

-1 

-1 

2 

2 

1 

-5 

2 

8 

9.87 

1 

-1 

0 

2 

0 

5 

0 

9 

14.83 

0 

0 

0 

2 

2 

-5 

2 

10 

29.93 

1 

-1 

0 

0 

1 

5 

-3 

11 

73.05 

0 

3 

2 

•2 

2 

-5 

2 

12 

187.66 

0 

0 

2 

-2 

3 

13 

-2 

13 

3230.13 

-1 

-1 

2 

-1 

2 

13 

-5 

14 

3231.50 

-1 

0 

1 

0 

1 

-15 

3 

15 

6164.10 

-1 

1 

0 

1 

1 

7 

-4 

16 

4.00 

3 

0 

2 

2 

2 

-1 

1 

17 

4.08 

1 

0 

2 

4 

2 

-2 

1 

18 

4.58 

4 

0 

2 

0 

2 

-3 

1 

19 

4.68 

2 

0 

2 

2 

1 

-2 

1 

20 

4.79 

0 

0 

2 

4 

1 

-1 

1 

21 

5.56 

1 

1 

2 

2 

2 

1 

-1 

22 

5.80 

-1 

0 

2 

4 

1 

-3 

1 

23 

5.90 

-1 

-1 

2 

4 

2 

-2 

1 

24 

6.73 

2 

1 

2 

0 

2 

4 

-2 

25 

6.82 

0 

0 

4 

0 

2 

2 

-1 

26 

6.85 

2 

0 

2 

0 

0 

3 

0 

27 

6.98 

1 

0 

2 

1 

2 

3 

-1 

28 

7.08 

0 

0 

2 

2 

0 

4 

0 

29 

7.13 

2 

0 

0 

2 

1 

-1 

1 

30 

7.23 

0 

-1 

2 

2 

1 

-4 

2 

31 

7.34 

.2 

0 

2 

4 

1 

-2 

1 

32 

7.39 

0 

0 

0 

4 

1 

-2 

1 

33 

8.68 

1 

0 

4 

-2 

2 

2 

-1 

34 

8.73 

3 

0 

2 

-2 

1 

2 

-1 

35 

8.90 

1 

1 

2 

0 

1 

4 

-2 

36 

9.05 

-1 

0 

4 

0 

1 

2 

-1 

37 

9.11 

0 

1 

2 

1 

2 

-2 

1 

38 

9.17 

-3 

0 

0 

0 

1 

-1 

1 

39 

9.33 

0 

0 

2 

1 

1 

3 

-1 

40 

9.35 

1 

-1 

2 

0 

1 

-4 

2 
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Table A. VI cont. 


ZMOA 1990-2 Theory of Nutation: Additional Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kj\ kj2 ^3 kj4 


•4o; Boj 

(0".00001) 

41 

9.60 

-1 

0 

0 

-2 

1 

-4 

3 

42 

10.07 

1 

0 

0 

-4 

1 

-1 

1 

43 

10.10 

-1 

0 

0 

4 

1 

-2 

1 

44 

10.37 

-1 

-1 

0 

4 

0 

1 

0 

45 

12.38 

2 

1 

2 

-2 

2 

3 

-1 

46 

12.64 

0 

0 

4 

-2 

1 

2 

-1 

47 

13.22 

1 

0 

2 

-1 

2 

-3 

1 

48 

13.28 

2 

1 

0 

0 

0 

-3 

0 

49 

13.69 

0 

0 

2 

0 

3 

2 

0 

50 

14.22 

0 

1 

0 

2 

1 

2 

-1 

51 

14.25 

1 

0 

0 

1 

0 

-3 

0 

52 

14.32 

2 

-1 

0 

0 

0 

4 

0 

53 

14.60 

-2 

0 

2 

2 

1 

2 

-1 

54 

14.70 

0 

0 

0 

-2 

2 

1 

-1 

55 

15.35 

0 

1 

0 

-2 

1 

-3 

2 

56 

15.42 

0 

-1 

0 

2 

1 

-2 

1 

57 

15.87 

2 

0 

0 

-4 

1 

-1 

1 

58 

15.94 

-2 

0 

0 

4 

1 

1 

-1 

59 

16.06 

0 

-2 

0 

2 

0 

2 

0 

60 

16.10 

0 

0 

2 

-4 

1 

-1 

1 

61 

22.40 

1 

1 

2 

-2 

1 

3 

-1 

62 

25.22 

-1 

1 

2 

0 

2 

4 

-2 

63 

25.53 

-1 

-1 

0 

0 

1 

2 

-1 

64 

25.72 

1 

1 

0 

0 

1 

-3 

2 

65 

26.77 

1 

0 

-2 

0 

1 

3 

-1 

66 

27.32 

0 

0 

1 

0 

1 

-2 

0 

67 

29.26 

-1 

-1 

2 

0 

2 

-2 

1 

68 

29.39 

-1 

1 

0 

2 

1 

-1 

1 

69 

29.40 

0 

0 

0 

-1 

1 

3 

-2 

70 

29.66 

0 

0 

0 

1 

1 

-4 

2 

71 

29.67 

-1 

1 

0 

0 

1 

-2 

2 

72 

31.52 

1 

0 

0 

-2 

2 

3 

-1 

73 

32.11 

-1 

0 

0 

2 

2 

-4 

2 

74 

32.45 

-1 

0 

2 

-2 

2 

3 

-1 

75 

35.80 

-1 

1 

2 

-2 

1 

-1 

1 

76 

38.52 

-1 

-2 

0 

2 

0 

3 

0 

77 

38.74 

1 

0 

2 

-4 

1 

-4 

2 

78 

121.75 

0 

3 

0 

0 

0 

3 

0 

79 

129.17 

-2 

-1 

0 

2 

1 

-2 

1 

80 

i 

177.85 

0 

-2 

0 

0 

1 

-1 

1 
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Table A. VI cont. 


ZMOA 1990-2 Theory of Nutation: Additional Terms 


Index 

j 

Period 

(days) 

Argument coefficient 
kji k }2 k } - 3 kj4 

/Cj5 

Aq] Bo j 

( 0 ". 00001 ) 

81 

219.17 

2 

0 

0 

-2 

2 

-3 

1 

82 

285.41 

-2 

1 

2 

0 

1 

-1 

0 

83 

297.91 

-2 

1 

2 

0 

2 

-1 

0 

84 

313.04 

-1 

0 

2 

-1 

1 

-4 

1 

85 

329.82 

0 

-1 

0 

0 

2 

4 

-1 

86 

438.33 

1 

0 

0 

.1 

1 

3 

-1 

87 

471.95 

-2 

-1 

2 

0 

2 

2 

-1 

88 

507.16 

-2 

-1 

2 

0 

1 

3 

0 

89 

552.62 

-3 

0 

2 

1 

2 

2 

-1 

90 

2266.13 

0 

0 

0 

0 

3 

-2 

0 

91 

6159.14 

-1 

0 

1 

0 

2 

3 

-1 


Table A. VII 

ZMOA 1990-2 Theory of Nutation: Out-of-Phase Terms 


Index 

j 

Period 

(days) 

kji 

Argument coefficient 
kj2 kj$ kj 4 

kjh 

A 2 ; B?) 

( 0 ". 00001 ) 

1 

6798.38 

0 

0 

0 

0 

1 

375 

168 

2 

182.62 

0 

0 

2 

-2 

2 

o 

in 

rH 

1 

-52 

3 

365.26 

0 

1 

0 

0 

0 

82 

-28 

4 

13.66 

0 

0 

2 

0 

2 

-2 

1 

5 

3399.19 

0 

0 

0 

0 

2 

-8 

-3 

6 

27.56 

1 

0 

0 

0 

0 

-3 

3 

7 

121.75 j 

0 

1 

2 

-2 

2 

-6 

-2 

8 

13.63 | 

0 

0 

2 

0 

1 

-1 

0 

9 

9.13 ! 

1 

0 

2 

0 

2 

3 

1 

10 

365.22 

1 

-1 

2 

-2 

2 

1 

1 

11 

31.81 

1 

0 

0 

-2 

0 

-1 

-1 

12 

177.84 

0 

0 

2 

-2 

1 

2 

0 

13 

27.09 

-1 

0 

2 

0 

2 

1 

0 

14 

386.00 

0 

1 

0 

0 

1 

1 

0 

15 

346.64 

0 

-1 

0 

0 

1 

1 

0 

16 

411.78 

1 

0 

0 

-1 

0 

-3 

1 
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Table A. VIII 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 
kjl kj2 kjZ k/4 

kjs 

A oy 
(Mas 

A n 

) 

Bo ) Bn 

(m qs ) 

1 

6798.38 

0 

0 

0 

0 

1 

-17199635 - 

17416 

9203279 

884 

2 

3399.19 

0 

0 

0 

0 

2 

207281 

21 

-89691 

47 

3 

1305.48 

-2 

0 

2 

0 

1 

4596 

5 

-2425 

-1 

4 

1095.18 

2 

0 

-2 

0 

0 

1100 

0 

11 

0 

5 

1615.75 

-2 

0 

2 

0 

2 

-312 

0 

131 

0 

6 

3232.86 

1 

-1 

0 

-1 

0 

-327 

0 

-1 

0 

7 

6786.32 

0 

-2 

2 

-2 

1 

-130 

0 

69 

0 

8 

943.23 

2 

0 

-2 

0 

1 

71 

0 

-37 

0 

9 

182.62 

0 

0 

2 

-2 

2 

-1318678 

-155 

573570 - 

305 

10 

365.26 

0 

1 

0 

0 

0 

142639 

-351 

5312 

-13 

11 

121,75 

0 

1 

2 

-2 

2 

-51735 

123 

22459 

-68 

12 

365.22 

0 

-1 

2 

-2 

2 

21686 

-50 

-9571 

30 

13 

177.84 

0 

0 

2 

-2 

1 

12851 

14 

-6900 

-1 

14 

205.89 

2 

0 

0 

-2 

0 

4768 

0 

43 

0 

15 

173.31 

0 

0 

2 

-2 

0 

-2177 

0 

-14 

0 

16 

182.63 

0 

2 

0 

0 

0 

1670 

-8 

13 

0 

17 

386.00 

0 

1 

0 

0 

1 

-1449 

-2 

841 

0 

18 

91.31 

0 

2 

2 

-2 

2 

-1579 

7 

685 

-4 

19 

346.64 

0 

-1 

0 

0 

1 

-1211 

1 

619 

0 

20 

199.84 

-2 

0 

0 

2 

1 

-576 

-1 

303 

0 

21 

346.60 

0 

-1 

2 

-2 

1 

-482 

-1 

270 

0 

22 

212.32 

2 

0 

0 

-2 

1 

402 

0 

-220 

0 

23 

119.61 

0 

1 

2 

-2 

1 

355 

0 

-194 

0 

24 

411.78 

1 

0 

0 

-1 

0 

-382 

0 

-28 

0 

25 

131.67 

2 

1 

0 

-2 

0 

115 

0 

1 

0 

26 

169.00 

0 

0 

-2 

2 

1 

87 

0 

-44 

0 

27 

329.79 

0 

1 

-2 

2 

0 

-88 

0 

2 

0 

28 

409.23 

0 

1 

0 

0 

2 

71 

0 

-32 

0 1 

29 

388.27 

-1 

0 

0 

1 

1 

88 

0 

-44 

0 

30 

117.54 

0 

1 

2 

-2 

0 

-66 

0 

0 

0 

31 

13.66 

0 

0 

2 

0 

2 

-227422 

-24 

97765 

-48 

32 

27.56 

1 

0 

0 

0 

0 

71162 

7 

-669 

0 

33 

13.63 

0 

0 

2 

0 

1 

-38700 

-37 

20065 

2 

34 

9.13 

1 

0 

2 

0 

2 

-30078 

-4 

12875 

-6 

35 

31.81 

1 

0 

0 

-2 

0 

-15710 

-1 

-123 

0 

36 

27.09 

-1 

0 

2 

0 

2 

12350 

1 

-5333 

3 

37 

14.77 

0 

0 

0 

2 

0 

6340 

1 

-121 

0 

38 

27.67 

1 

0 

0 

0 

1 

6311 

6 

-3323 

0 

39 

27.44 

-1 

0 

0 

0 

1 

-5807 

-6 

3144 

0 

40 

9.56 

-1 

0 

2 

2 

2 

-5951 

-1 

2550 

-1 
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Table A. VIII cont. 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 

fcyi k,2 k, 3 k }< 

kjb 

Aoy 

(Mas) 


Bo, 

(Mas 

Bi, 

41 

9.12 

1 

0 

2 

0 

1 

-5147 

-4 

2631 

0 

42 

7.10 

0 

0 

2 

2 

2 

-3842 

0 

1638 

-1 

43 

13.78 

2 

0 

0 

0 

0 

2925 

0 

-60 

0 

44 

23.94 

1 

0 

2 

-2 

2 

2861 

0 

-1234 

1 

45 

6.86 

2 

0 

2 

0 

2 

-3091 

0 

1319 

-1 

46 

13.61 

0 

0 

2 

0 

0 

2586 

0 

-54 

0 

47 

26.98 

-1 

0 

2 

0 

1 

2042 

2 

-1076 

0 

48 

31.96 

-1 

0 

0 

2 

1 

1516 

1 

-801 

0 

49 

31.66 

1 

0 

0 

-2 

1 

-1288 

-1 

695 

0 

50 

9.54 

-1 

0 

2 

2 

1 

-1018 

-1 

521 

0 

51 

34.85 

1 

1 

0 

-2 

0 

-735 

0 

-5 

0 

52 

13.17 

0 

1 

2 

0 

2 

756 

-2 

-325 

0 

53 

14.19 

0 

-1 

2 

0 

2 

-713 

2 

307 

0 

54 

5.64 

1 

0 

2 

2 

2 

-765 

0 

324 

0 

55 

9.61 

1 

0 

0 

2 

0 

658 

0 

-20 

0 

56 

12.81 

2 

0 

2 

-2 

2 

642 

0 

-277 

0 

57 

14.80 

0 

0 

0 

2 

1 

-630 

-1 

327 

0 

58 

7.09 

0 

0 

2 

2 

1 

-661 

-1 

334 

0 

59 

23.86 

1 

0 

2 

-2 

1 

571 

1 

-307 

0 

60 

14.73 

0 

0 

0 

-2 

1 

-494 

-1 

273 

0 

61 

29.80 

1 

-1 

0 

0 

0 

471 

0 

-4 

0 

62 

6.85 

2 

0 

2 

0 

1 

-531 

o 

268 

0 

63 

15.39 

0 

1 

0 

-2 

0 

-435 

o 

-8 

0 

64 

26.88 

1 

0 

-2 

0 

0 

404 

0 

4 

0 

65 

29.53 

0 

0 

0 

1 

0 

-422 

0 

3 

0 

66 

25.62 

1 

1 

0 

0 

0 

-338 

0 

3 

0 

67 

9.11 

1 

0 

2 

0 

0 

334 

0 

-11 

0 

68 

9.37 

1 

-1 

2 

0 

2 

-287 

0 

■ 

0 

69 

9.81 

-1 

-1 

2 

2 

2 

-281 

0 


0 

70 

13.75 

-2 

0 

0 

0 

1 

-229 

0 

126 

0 

71 

5.49 

3 

0 

2 

0 

2 

-288 

0 


0 

72 

7.24 

0 

-1 

2 

2 

2 

-263 

0 

114 

0 

73 

8.91 

1 

1 

2 

0 

2 

246 

0 

-106 

0 

74 

32.61 

-1 

0 

2 

-2 

1 

-200 

0 

107 

0 

75 

13.81 

2 

0 

0 

0 

1 

218 

0 

-114 

0 

76 

27.78 

1 

0 

0 

0 

2 

-197 

0 

85 

0 

77 

9.18 

3 

0 

0 

0 

0 

157 

0 

-4 

0 

78 

9.34 

0 

0 

2 

1 

2 

165 

0 

-72 

0 

79 

27.33 

-1 

0 

0 

0 

2 

140 

0 

-61 

0 

80 

10.08 

1 

0 

0 

-4 

0 

-133 

0 

-4 

0 

i 
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Table A. VIII cont. 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 
fcyi ^;3 ^j4 ^j5 

Aoj . 

(jxas) 


Bqj 

(Mas 

Bn 

) 

81 

14.63 

-2 

0 

2 

2 

2 

139 

0 

-60 

0 

82 

5.80 

-1 

0 

2 

4 

2 

.151 

0 

65 

0 

83 

15.91 

2 

0 

0 

-4 

0 

-128 

0 

1 

0 

84 

22.47 

1 

1 

2 

-2 

2 

128 

0 

-55 

0 

85 

5.64 

1 

0 

2 

2 

1 

-132 

0 

66 

0 

86 

7.35 

-2 

0 

2 

4 

2 

-121 

0 

52 

0 

87 

9.06 

-1 

0 

4 

0 

2 

115 

0 

-49 

0 

88 

29.26 

1 

-1 

0 

-2 

0 

95 

0 

0 

0 

89 

12.79 

2 

0 

2 

-2 

1 

101 

0 

-54 

0 

90 

4.68 

2 

0 

2 

2 

2 

-107 

0 

47 

0 

91 

9.63 

1 

0 

0 

2 

1 

-95 

0 

49 

0 

92 

12.66 

0 

0 

4 

-2 

2 

91 

0 

-39 

0 

93 

8.75 

3 

0 

2 

-2 

2 

93 

0 

-40 

0 

94 

23.77 

1 

0 

2 

-2 

0 

-77 

0 

0 

0 

95 

13.14 

0 

1 

2 

0 

1 

81 

0 

-42 

0 

96 

35.03 

-1 

-1 

0 

2 

1 

75 

0 

-39 

0 

97 

13.58 

0 

0 

-2 

0 

1 

-68 

0 

36 

0 

98 

25.42 

0 

0 

2 

-1 

2 

-67 

0 

31 

0 

99 

14.19 

0 

1 

0 

2 

0 

-67 

0 

-1 

0 

100 

9.53 

1 

0 

-2 

-2 

0 

-64 

0 

1 

0 

101 

14.16 

0 

-1 

2 

0 

1 

-65 

0 

35 

0 

102 

34.67 

1 

1 

0 

-2 

1 

-61 

0 

32 

0 

103 

32.76 

1 

0 

-2 

2 

0 

-58 

0 

0 

0 

104 

7.13 

2 

0 

0 

2 

0 

59 

0 

1 

0 

105 

4.79 

0 

0 

2 

4 

2 

-68 

0 

29 

0 

106 

27.32 

0 

1 

0 

1 

0 

57 

0 

0 

0 

107 

66068.15 

-1 

1 

0 

1 

2 

-52 

0 

23 

0 

108 

65502.23 

-1 

0 

1 

0 

3 

-6 

0 

0 

0 

109 

6164.10 

-1 

1 

0 

1 

1 

126 

0 

-68 

0 

110 

6159.14 

-1 

0 

1 

0 

2 

36 

0 

18 

0 

111 

3396.17 

0 

2 

-2 

2 

0 

-9 

0 

0 

0 

112 

3231.50 

-1 

0 

1 

0 

1 

-103 

0 

-89 

0 

113 

2266.13 

0 

0 

0 

0 

3 

-22 

0 

10 

0 

114 

2190.35 

-1 

0 

1 

0 

0 

33 

0 

0 

0 

115 

2189.72 

1 

1 

-2 

I 

-1 

-7 

0 

0 

0 

116 

1616.43 

-2 

2 

0 

2 

0 

-22 

0 

0 

0 

117 

548.04 

-2 

-1 

2 

0 

0 

-4 

0 

0 

0 

118 

507.16 

-2 

-1 

2 

0 

1 

6 

0 

-4 

0 

119 

471.89 

-2 

1 

0 

2 

0 

1 

0 

4 

0 

120 

438.33 

-1 

0 

0 

1 

-1 

-9 

0 

-8 

0 
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Table A. VIII cont. 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 



Period 

(days) 

k } - 1 

Argument coefficient 

^j2 ^j3 


Aqj 

(pas) 

Ms 

•Boy 

(/ias 

Bis 

) 

121 

329.82 

0 

-1 

0 

0 

2 

57 

0 

-25 

0 

122 

313.04 

1 

0 

-2 

1 

-1 

23 

0 

14 

0 

123 

299.26 

1 

0 

-2 

1 

0 

-7 

0 

0 

0 

124 

285.41 

-2 

1 

2 

0 

1 

9 

0 

0 

0 

125 

219.17 

-2 

0 

0 

2 

-2 

11 

0 

0 

0 

126 

194.13 

-2 

0 

0 

2 

2 

17 

0 

-8 

0 

127 

187.67 

0 

2 

0 

0 

1 

-9 

0 

5 

0 

128 

187.66 

0 

0 

-2 

2 

-3 

-122 

0 

-28 

0 

129 

177.85 

0 

-2 

0 

0 

1 

-6 

0 

0 

0 

130 

134.27 

-2 

-1 

0 

2 

-1 

-8 

0 

0 

0 

131 

129.17 

-2 

-1 

0 

2 

1 

-13 

0 

7 

0 

132 

126.51 

1 

0 

2 

-3 

2 

6 

0 

0 

0 

133 

95.42 

-2 

0 

-2 

4 

-1 

-10 

0 

-5 

0 

134 

90.10 

0 

-2 

-2 

2 

-1 

-8 

0 

-5 

0 

135 

73.05 

0 

3 

2 

-2 

2 

-44 

0 

19 

0 

136 

38.96 

-1 

0 

-2 

4 

0 

-7 

0 

0 

0 

137 

38.74 

-1 

0 

-2 

4 

-1 

47 

0 

24 

0 

138 

38.52 

-1 

-2 

0 

2 

0 

27 

0 

0 

0 

139 

38.52 

-1 

0 

-2 

4 

-2 

-7 

0 

0 

0 

140 

37.63 

-3 

0 

0 

4 

0 

16 

o 

0 

0 

141 

35.80 

1 

-1 

-2 

2 

-1 

11 

0 

5 

0 

142 

35.23 

0 

0 

2 

-3 

2 

-20 

0 

8 

0 

143 

34.48 

-2 

0 

0 

3 

0 

-5 

0 

0 

0 

144 

32.45 

1 

0 

-2 

2 

-2 

-36 

0 

-16 

0 

145 

32.45 

1 

-2 

0 

0 

0 

8 

0 

0 

0 

146 

32.11 

-1 

0 

0 

2 

2 

-38 

0 

17 

0 

147 

31.52 

-1 

0 

0 

2 

-2 

-32 

0 

-14 

0 

148 

29.93 

1 

-1 

0 

0 

1 

51 

0 

-27 

0 

149 

29.67 

1 

-1 

0 

0 

-1 

31 

0 

17 

0 

150 

29.66 

0 

0 

0 

1 

1 

-37 

0 

20 

0 

151 

29.40 

0 

0 

0 

1 

-1 

-33 

0 

-18 

0 

152 

29.39 

-1 

1 

0 

2 

1 

-10 

0 

6 

0 

153 

29.26 

-1 

-1 

2 

0 

2 

-16 

0 

6 

0 

154 

29.14 

-1 

1 

0 

2 

-1 

-6 

0 

0 

0 

155 

29.14 

-1 

-1 

2 

0 

1 

-6 

0 

0 

0 

156 

28.15 

-3 

0 

2 

0 

1 

-8 

0 

0 

0 

157 

27.43 

0 

1 

0 

1 

1 

5 

0 

0 

0 

158 

27.32 

0 

1 

-2 

1 

-2 

-9 

0 

0 

0 

159 

27.32 

0 

0 

1 

0 

1 

-17 

0 

-15 

0 

160 

27.21 

0 

0 

1 

0 

0 

8 

0 

0 

0 
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Table A. VIII cont. 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 


Index 

j 

Period 

(days) 


Argument coefficient 

k j2 k , 3 k , 4 

£ j 5 

•Aoy 

(Mas) 


B oj 
(/uas 

B iy 

161 

27.09 

-1 

2 

0 

2 

0 

-8 

0 

0 

0 

162 

26.77 

-1 

0 

2 

0 

-1 

-27 

0 

-14 

0 

163 

25.72 

1 

1 

0 

0 

1 

-34 

0 

19 

0 

164 

25.62 

-1 

1 

-2 

2 

-2 

7 

0 

0 

0 

165 

25.53 

1 

1 

0 

0 

-1 

-21 

0 

-12 

0 

166 

25.33 

0 

0 

-2 

1 

-1 

10 

0 

5 

0 

167 

25.22 

-1 

1 

2 

0 

2 

45 

0 

-20 

0 

168 

25.13 

-1 

1 

2 

0 

1 

9 

0 

-5 

0 

169 

23.43 

1 

0 

-4 

2 

-2 

9 

0 

0 

0 

170 

22.40 

-1 

-1 

-2 

2 

-1 

-23 

0 

-12 

0 

171 

16.63 

-2 

-1 

0 

4 

0 

12 

0 

0 

0 

172 

16.10 

0 

0 

-2 

4 

-1 

16 

0 

9 

0 

173 

16.06 

0 

-2 

0 

2 

0 

21 

0 

0 

0 

174 

15.94 

-2 

0 

0 

4 

1 

11 

0 

-6 

0 

175 

15.87 

-2 

0 

0 

4 

-1 

11 

0 

5 

0 

176 

15.42 

0 

-1 

0 

2 

1 

-14 

0 

7 

0 

177 

15.35 

0 

-1 

0 

2 

-1 

40 

0 

18 

0 

178 

15.31 

-1 

0 

0 

3 

0 

-11 

0 

0 

0 

179 

15.24 

-2 

-1 

2 

2 

2 

6 

0 

0 

0 

180 

14.93 

2 

0 

-2 

2 

-1 

8 

0 

0 

0 

181 

14.83 

0 

0 

0 

2 

2 

-46 

0 

20 

0 

182 

14.77 

0 

-2 

2 

0 

2 

-8 

0 

0 

0 

183 

14.70 

0 

0 

0 

2 

-2 

-13 

0 

-5 

0 

184 

14.60 

-2 

0 

2 

2 

1 

21 

0 

-11 

0 

185 

14.57 

-2 

0 

2 

2 

0 

-6 

0 

0 

0 

186 

14.32 

2 

-1 

0 

0 

0 

38 

0 

0 

0 

187 

14.25 

1 

0 

0 

1 

0 

-33 

0 

0 

0 

188 

14.22 

0 

1 

0 

2 

1 

18 

0 

-9 

0 

189 

14.16 

0 

1 

0 

2 

-1 

-5 

0 

0 

0 

190 

14.13 

-1 

0 

2 

1 

2 

8 

0 

0 

0 

191 

13.83 

2 

0 

0 

0 

2 

-8 

0 

0 

0 

192 

13.72 

2 

0 

0 

0 

-2 

-6 

0 

0 

0 

193 

13.72 

1 

1 

0 

1 

0 

5 

0 

0 

0 

194 

13.69 

0 

0 

2 

0 

3 

20 

0 

0 

0 

195 

13.49 

-2 

0 

4 

0 

2 

-7 

0 

0 

0 

196 

13.28 

2 

1 

0 

0 

0 

-28 

0 

0 

0 

197 

13.28 

-2 

1 

-2 

2 

-2 

-5 

0 

0 

0 

198 

13.22 

-1 

0 

-2 

1 

-2 

33 

0 

13 

0 

199 

12.76 

-2 

0 

-2 

2 

0 

8 

0 

0 

0 

200 

12.71 

0 

2 

2 

0 

2 

7 

0 

0 

0 
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Table A. VIII cont. 


Kinoshita-Souchay Non-rigid Earth Theory of Nutation 



Index 

j 

Period 

(days) 

k 3 1 

Argument coefficient 

^j'2 

*ys 

-^Oj Aij 

( fias ) 

B 0] Bij 

(Mas) 

201 

12.64 

0 

0 

-4 

2 

-1 

-12 

0 

-7 

0 

202 

12.38 

-2 

-1 

-2 

2 

-2 

-23 

0 

-10 

0 

203 

10.37 

-1 

-1 

0 

4 

0 

15 

0 

0 

0 

204 

10.10 

-1 

0 

0 

4 

1 

-18 

0 

9 

0 

205 

10.08 

-1 

-2 

2 

2 

2 

-9 

0 

0 

0 

206 

10.07 

-1 

0 

0 

4 

-1 

12 

0 

6 

0 

207 

9.87 

1 

-1 

0 

2 

0 

48 

0 

1 

0 

208 

9.80 

-1 

-1 

2 

2 

1 

-46 

0 

24 

0 

209 

9.60 

1 

0 

0 

2 

-1 

50 

0 

27 

0 

210 

9.37 

1 

1 

0 

2 

0 

-11 

0 

0 

0 

211 

9.35 

1 

-1 

2 

0 

1 

-41 

0 

22 

0 

212 

9.33 

0 

0 

2 

1 

1 

27 

0 

-14 

0 

213 

9.31 

-1 

1 

2 

2 

2 

56 

0 

-23 

0 

214 

9.30 

-1 

1 

2 

2 

1 

7 

0 

0 

0 

215 

9.20 

3 

0 

0 

0 

1 

7 

0 

0 

0 

216 

9.17 

3 

0 

0 

0 

-1 

14 

0 

8 

0 

217 

9.11 

0 

1 

2 

1 

2 

-24 

0 

10 

0 

218 

9.10 

1 

0 

2 

0 

-1 

10 

0 

5 

0 

219 

9.05 

*1 

0 

4 

0 

1 

19 

0 

-9 

0 

220 

8.93 

‘ 2 

0 

-2 

1 

-2 

6 

0 

0 

0 

221 

8.90 

1 

1 

2 

0 

1 

36 

0 

-18 

0 

222 

8.73 

-3 

0 

-2 

2 

-1 

-12 

0 

-7 

0 

223 

8.68 

-1 

0 

•4 

2 

-2 

-19 

0 

-9 

o 

224 

7.53 

0 

-1 

0 

4 

0 

6 

0 

0 

0 

225 

7.50 

-2 

-1 

2 

4 

2 

-9 

0 

0 

o 

226 

7.39 

0 

0 

0 

4 

1 

-13 

0 

7 

0 

227 

7.38 

0 

0 

0 

4 

0 

48 

0 

1 

o 

228 

7.38 

0 

-2 

2 

2 

2 

-11 

0 

5 

0 

229 

7.34 

-2 

0 

2 

4 

1 

-14 

0 

10 

o 

230 

7.23 

0 

-1 

2 

2 

1 

-44 

0 

23 

0 

231 

7.22 

-1 

0 

2 

3 

2 

8 

0 

0 

o 

232 

7.13 

2 

0 

0 

2 

1 

-9 

0 

5 

0 

233 

7.12 

2 

0 

0 

2 

-1 

5 

0 

0 

o 

234 

7.08 

0 

0 

2 

2 

0 

41 

0 

1 

o 

235 

6.99 

2 

-1 

2 

0 

2 

-48 

0 

21 

0 

236 

6.98 

1 

0 

2 

1 

2 

33 

0 

-14 

0 

237 

6.98 

2 

-1 

2 

0 

1 

-6 

0 

0 

o 

238 

6.97 

1 

0 

2 

1 

1 

5 

0 

0 

o 

239 

6,96 

0 

1 

2 

2 

2 

54 

0 

-16 

o 

240 

6.95 

i 

0 

1 

2 

2 

1 

7 

0 

0 

0 
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Kinoshita-Souchay Non-rigid Earth Theory of Nutation: Planetary Terms 
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Table A.X 


Woolard Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 
kji kj 2 kj 4 


•Aoy 

(0".0001) 

Boy Bij 
( 0 ". 0001 ) 

1 

6798.4 

0 

0 

0 

0 

1 

-172327 

-173.7 

92100 

9.1 

2 

3399.2 

0 

0 

0 

0 

2 

2088 

0.2 

-904 

0.4 

3 

1305.5 

-2 

0 

2 

0 

1 

45 

0.0 

-24 

0.0 

4 

1095.2 

2 

0 

-2 

0 

0 

10 

0.0 

0 

0.0 

5 

1615.7 

-2 

0 

2 

0 

2 

-3 

0.0 

2 

0.0 

6 

3232.9 

1 

-1 

0 

-1 

0 

-2 

0.0 

0 

0.0 

7 

6786.3 

0 

-2 

2 

-2 

1 

-4 

0.0 

2 

0.0 

8 

182.6 

0 

0 

2 

-2 

2 

-12729 

-1.3 

5522 

-2.9 

9 

365.3 

0 

1 

0 

0 

0 

1261 

-3.1 

0 

0.0 

10 

121.7 

0 

1 

2 

-2 

2 

-497 

1.2 

216 

- 0.6 

11 

365.2 

0 

-1 

2 

-2 

2 

214 

-0.5 

-93 

0.3 

12 

177.8 

0 

0 

2 

-2 

1 

124 

0.1 

-66 

0.0 

13 

205.9 

2 

0 

0 

-2 

0 

45 

0.0 

0 

0.0 

14 

173.3 

0 

0 

2 

-2 

0 

-21 

0.0 

0 

0.0 

15 

182.6 

0 

2 

0 

0 

0 

16 

- 0.1 

0 

0.0 

16 

386.0 

0 

1 

0 

0 

1 

-15 

0.0 

8 

0.0 

17 

91.3 

0 

2 

2 

-2 

2 

-15 

0.1 

7 

0.0 

18 

346.6 

0 

-1 

0 

0 

1 

-10 

0.0 

5 

0.0 

19 

199.8 

-2 

0 

0 

2 

1 

-5 

0.0 

3 

0.0 

20 

346.6 

0 

-1 

2 

-2 

1 

-5 

0.0 

3 

0.0 

21 

212.3 

2 

0 

0 

-2 

1 

4 

0.0 

-2 

0.0 

22 

119.6 

0 

1 

2 

-2 

1 

3 

0.0 

-2 

0.0 

23 

411.8 

1 

0 

0 

-1 

0 

-3 

0.0 

0 

0.0 

24 

13.7 

0 

0 

2 

0 

2 

-2037 

- 0.2 

884 

-0.5 

25 

27.6 

1 

0 

0 

0 

0 

675 

0.1 

0 

0.0 

26 

13.6 

0 

0 

2 

0 

1 

-342 

-0.4 

183 

0.0 

27 

9.1 

1 

0 

2 

0 

2 

-261 

0.0 

113 

- 0.1 

28 

31.8 

1 

0 

0 

-2 

0 

-149 

0.0 

0 

0.0 

29 

27.1 

-1 

0 

2 

0 

2 

114 

0.0 

-50 

0.0 

30 

14.8 

0 

0 

0 

2 

0 

60 

0.0 

0 

0.0 

31 

27.7 

1 

0 

0 

0 

1 

58 

0.0 

-31 

0.0 

32 

27.4 

-1 

0 

0 

0 

1 

-57 

0.0 

30 

0.0 

33 

9.6 

-1 

0 

2 

2 

2 

-52 

0.0 

22 

0.0 

34 

9.1 

1 

0 

2 

0 

1 

-44 

0.0 

23 

0.0 

35 

7.1 

0 

0 

2 

2 

2 

-32 

0.0 

14 

0.0 
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Table A.X cont. 


Woolard Theory of Nutation 


Index 

j 

Period 

(days) 

Argument coefficient 

^'3 



B 0 j Bij 

(0".0001) 

36 

13.8 

2 

0 

0 

0 

0 

28 

0.0 

0 

0.0 

37 

23.9 

1 

0 

2 

-2 

2 

26 

0.0 

-11 

0.0 

38 

6.9 

2 

0 

2 

0 

2 

-26 

0.0 

11 

0.0 

39 

13.6 

0 

0 

2 

0 

0 

25 

0.0 

0 

0.0 

40 

27.0 

-1 

0 

2 

0 

1 

19 

0.0 

-10 

0.0 

41 


-1 

0 

0 

2 

1 

14 

0.0 

-7 

0,0 

42 

31.7 

1 

0 

0 

-2 

1 

-13 

0.0 

7 

0.0 

43 

9.5 

-1 

0 

2 

2 

1 

-9 

0.0 

5 

0.0 

44 

34.8 

1 

1 

0 

-2 

0 

-7 

0.0 

0 

0.0 

45 

13.2 

0 

1 

2 

0 

2 

7 

0.0 

-3 

0.0 

46 

14.2 

0 

-1 

2 

0 

2 

-6 

0.0 

3 

0.0 

47 

5.6 

1 

0 

2 

2 

2 

-6 

0.0 

3 

0.0 

48 

9.6 

1 

0 

0 

2 

0 

6 

0.0 

0 

0.0 

49 

12.8 

2 

0 

2 

-2 

2 

6 

0.0 

-2 

0.0 


14.8 

0 

0 

0 

2 

1 

-6 

0.0 

3 

0.0 

51 

7.1 

0 

0 

2 

2 

1 

-5 

0.0 

3 

0.0 

52 

23.9 

i 

0 

2 

-2 

1 

5 

0.0 

-3 

0.0 

53 

14.7 

0 

0 

0 

-2 

1 

-5 

0.0 

3 

0.0 

54 

29.8 

1 

-1 

0 

0 

0 

4 

0.0 

0 

0.0 

55 

6.9 

2 

0 

2 

0 

1 

-4 

0.0 

2 

0.0 

56 

15.4 

0 

1 

0 

-2 

0 

-4 

0.0 

0 

0.0 

57 

26.9 

1 

0 

-2 

0 

0 

4 

0.0 

0 

0.0 

58 

29.5 

0 


0 

1 

0 

-4 

0.0 

0 

0.0 

59 

25.6 

1 

1 

0 

0 

0 

-3 

0.0 

0 

0.0 

60 

9.1 

1 

0 

2 

0 

0 

3 

0.0 

0 

0.0 

■■ 

9.4 

1 

-1 

2 

0 

2 

-3 

0.0 

0 

0.0 

62 

9.8 

.1 

-1 

2 

2 

2 

-2 

0.0 

0 

0.0 

63 

13.7 

-2 

0 

0 

0 

1 

-2 

0.0 

0 

0.0 

64 

5.5 

3 

0 

2 

0 

2 

-2 

0.0 

0 

0.0 

65 

7.2 

0 

-1 

2 

2 

2 

-2 

0.0 

0 

0.0 

66 

8.9 

1 

1 

2 

0 

2 

2 

0.0 

0 

0.0 

67 

32.6 

-1 

0 

2 

-2 

1 

-2 

0.0 

0 

0.0 

68 

13.8 

2 

0 

0 

0 

1 

2 

0.0 

0 

0.0 

69 

27.8 

1 

0 

0 

0 

2 

-2 

0.0 

0 

0.0 
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APPENDIX B 

GLOSSARY OP “MODEST” PARAMETERS 

For the convenience of users of MODEST, Table B.I identifies the names of adjustable parameters in 
the code with the notation of this document. Brief definitions and either references to equations (in 
parentheses) or sections (no parentheses) are also given. 


Table B.I 

Glossary of MODEST Parameters 


Parameter 

MODEST name 

Definition 

Reference 


RSPINAX aaaaaaaa 

Cylindrical 

(2.36) 

A 

LDNGTUD aaaaaaaa 

station 

(2.37) 

z 

P0LPR0J aaaaaaaa 

coordinates 

(2.38) 



Time rates of 

(2.36) 

A 

DLQN/DT aaaaaaaa 

change of 

(2.37) 

z 


station coordinates 

(2.38) 

X 

X aaaaaaaa 

Cartesian 

(2.39) 

y 

Y aaaaaaaa 

station 

(2.40) 

z 

Z aaaaaaaa 

coordinates 

(2.41) 

X 

DX/DT aaaaaaaa 

Time rates of 

(2.39) 

y 


change of 

(2.40) 

z 


station coordinates 

(2.41) 

l 

AXISOFF aaaaaaaa 

Antenna offset 

(2.199) 

h, l 

*L0VE # aaaaaaaa 

Love numbers 

(2.49) to (2.51) 

v, 

TIDEPHZ aaaaaaaa 

Tide lag 

(2.46) 


OLAMgCcc aaaaaaaa 
OLAMgScc aaaaaaaa 

Ocean loading 
amplitudes 

(2.74) 

f 

ATMOLOADaaaaaaaa 

Loading factor 

(2.80) 

"1rr N 

GEN REL GAMMA FACTOR 

PPN gamma 

(2.16) 

Cx 

RIGHT ASCEN . S8S88BSBB8B8 

Source RA 

(2.227) 

s 

DECLINATION ssssssssssss 

Source declination 

(2.227) 

a 

DRASCEN/DT ssssssssssss 

Time rates of change 

(2.91) 

6 

DDECLIN/DT ssssssssssss 

of RA, declination 

(2.92) 

K 

FLUX RATIO ssssssssssss 

Source 

(2.93) 

s 

SEPARATION bsssbsssbsss 

structure 

to 

e 

ORIENTATION ssssssssssss 

parameters 

(2.97) 


aaaaaaaa station name 

ssssssssssss source name 

cc component name (alphanumeric) 

g geometric direction (U, N, W) 

* V or H 
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Table B.I cont. 

Glossary of MODEST Parameters 


Parameter 

MODEST name 

Definition 

Reference 

©1,2 

t POLE MOTION 

Pole position 

(2.101-102) 

UT 1 - UTC 

UT1 MINUS UTC 

UTl - UTC 

2.6.1 

©1,2 

f POLE RATE 

Time rates of 

(2.111) 

©3 

UT1-UTC RATE 

change of UTPM 



UT1 TIDAL AMPLITUDE cnnn 

Tidal terms 

(2.114) 


tPOL TIDAL AMPLITUDE cnnn 

in UTPM 


6&x,v,e 

t AXIS TWEAK OFFSET 

Perturbation 

(2.166) 

$Qx,u,e 

t AXIS TWEAK RATE 

coefficients 

(2.166) 

Pls 

LUNI-SOLAR PRECESSION 

Precession 

(2.156) 

PPL 

PLANETARY PRECESSION 

constants 

(2.156) 

A 0 j 

NUTATION AMPLTD PSI cjjj 

Nutation 

(2.128) to 

Ai } 

NUTATION AMPLTD PSITcjjj 

amplitudes 

(2.133) 

Aij.zj 

NUTATION AMPLTD PSIA 



Bo, 

NUTATION AMPLTD EPS cjjj 



Bis 

NUTATION AMPLTD EPSTcjjj 



B 2S,3] 

NUTATION AMPLTD EPSA 




t X or Y 
t X, Y, or Z 

c component: S, C for sine, cosine 

j j j nutation series term number 

nnn tidal series term number 
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Table B.I cont. 

Glossary of MODEST Parameters 


Parameter 

MODEST name 

Definition 

Reference 

T c i 


Coefficients 

(3.1) 

Tc2 

C RATE aaaaaaaa 

in clock 

(3.1) 

T C 3 

DCRAT/DTaaaaaaaa 

model for 

(3.1) 

T C 4 

F OFFSETaaaaaaaa 

delay and 

(3.2) 

T c 5 

F DRIFT aaaaaaaa 

delay rate 

(3.2) 

Tc6 

C EPOS/Xaaaaaaaa 

S/X clock offset 

(3.3) 


DRYZTROPaaaaaaaa 

Dry zenith delay 

(4.3) 

PZwt 
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